2012-02-13 132 views
3

我正在R(glm)中运行逻辑回归。然后我设法绘制结果。我的代码如下:绘制R逻辑回归中的两条曲线

temperature.glm = glm(Response~Temperature, data=mydata,family=binomial) 

plot(mydata$Temperature,mydata$Response, ,xlab="Temperature",ylab="Probability of Response") 
curve(predict(temperature.glm,data.frame(Temperature=x),type="resp"),add=TRUE, col="red") 
points(mydata$Temperature,fitted(temperature.glm),pch=20) 
title(main="Response-Temperature with Fitted GLM Logistic Regression Line") 

我的问题是:

  1. 我怎么能在一个图形绘制2条Logistic回归曲线?
  2. 我从其他统计软件中得到了这两个系数。我怎样才能创建随机数据,插入这两套coef(Set 1和Set 2),然后产生两个逻辑回归曲线?

该机型:

    SET 1 
(Intercept)  -88.4505 
Temperature  2.9677 

        SET 2 
(Intercept) -88.585533 
Temperature  2.972168 

mydata是在2列和〜700行。

Response Temperature 
1 29.33 
1 30.37 
1 29.52 
1 29.66 
1 29.57 
1 30.04 
1 30.58 
1 30.41 
1 29.61 
1 30.51 
1 30.91 
1 30.74 
1 29.91 
1 29.99 
1 29.99 
1 29.99 
1 29.99 
1 29.99 
1 29.99 
1 30.71 
0 29.56 
0 29.56 
0 29.56 
0 29.56 
0 29.56 
0 29.57 
0 29.51 
+0

为什么你就不能调用'curve'(或'lines')的两倍,与不同的曲线值? – 2012-02-13 10:59:20

+0

另外,如果您提供可重现的数据集,回答您的问题会更容易。在这种情况下,我们无法访问'mydata',这会让事情变得更加困难。 – 2012-02-13 11:00:47

+0

最后,删除了你的签名。如果你想让人们知道你是埃迪,请将你的名字写在你的个人资料中。欢迎来到SO,顺便说一句。 – 2012-02-13 11:02:50

回答

15
  1. 要绘制一条曲线,你只需要定义响应和预测之间的关系,并指定要为其该曲线绘制的预测值的范围。例如: -

    dat <- structure(list(Response = c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
        1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 0L, 0L, 0L, 0L, 0L, 
        0L, 0L), Temperature = c(29.33, 30.37, 29.52, 29.66, 29.57, 30.04, 
        30.58, 30.41, 29.61, 30.51, 30.91, 30.74, 29.91, 29.99, 29.99, 
        29.99, 29.99, 29.99, 29.99, 30.71, 29.56, 29.56, 29.56, 29.56, 
        29.56, 29.57, 29.51)), .Names = c("Response", "Temperature"), 
        class = "data.frame", row.names = c(NA, -27L)) 
    
    temperature.glm <- glm(Response ~ Temperature, data=dat, family=binomial) 
    
    plot(dat$Temperature, dat$Response, xlab="Temperature", 
        ylab="Probability of Response") 
    curve(predict(temperature.glm, data.frame(Temperature=x), type="resp"), 
         add=TRUE, col="red") 
    # To add an additional curve, e.g. that which corresponds to 'Set 1': 
    curve(plogis(-88.4505 + 2.9677*x), min(dat$Temperature), 
         max(dat$Temperature), add=TRUE, lwd=2, lty=3) 
    legend('bottomright', c('temp.glm', 'Set 1'), lty=c(1, 3), 
         col=2:1, lwd=1:2, bty='n', cex=0.8) 
    

    在上面的第二curve电话,我们都在讲逻辑函数定义xy之间的关系。 plogis(z)的结果与评估1/(1+exp(-z))时的结果相同。参数min(dat$Temperature)max(dat$Temperature)定义了x的范围,其中y应该被评估。我们不需要告诉x指的是温度;当我们指定应该为预测值的范围评估响应时,这是隐含的。

    Adding additional curves to a plot

  2. 正如可以看到的,curve功能可以绘制的曲线,而无需模拟预测器(例如,温度)的数据。如果你仍然需要这样做,例如绘制符合特定型号的伯努利试验的一些模拟的结果,那么可以尝试以下方法:

    n <- 100 # size of random sample 
    
    # generate random temperature data (n draws, uniform b/w 27 and 33) 
    temp <- runif(n, 27, 33) 
    
    # Define a function to perform a Bernoulli trial for each value of temp, 
    # with probability of success for each trial determined by the logistic 
    # model with intercept = alpha and coef for temperature = beta. 
    # The function also plots the outcomes of these Bernoulli trials against the 
    # random temp data, and overlays the curve that corresponds to the model 
    # used to simulate the response data. 
    sim.response <- function(alpha, beta) { 
        y <- sapply(temp, function(x) rbinom(1, 1, plogis(alpha + beta*x))) 
        plot(y ~ temp, pch=20, xlab='Temperature', ylab='Response') 
        curve(plogis(alpha + beta*x), min(temp), max(temp), add=TRUE, lwd=2)  
        return(y) 
    } 
    

    实例:

    # Simulate response data for your model 'Set 1' 
    y <- sim.response(-88.4505, 2.9677) 
    
    # Simulate response data for your model 'Set 2' 
    y <- sim.response(-88.585533, 2.972168) 
    
    # Simulate response data for your model temperature.glm 
    # Here, coef(temperature.glm)[1] and coef(temperature.glm)[2] refer to 
    # the intercept and slope, respectively 
    y <- sim.response(coef(temperature.glm)[1], coef(temperature.glm)[2]) 
    

    下图显示了由上面的第一实施例制造的积即对温度的随机向量的每个值进行单个伯努利试验的结果,以及描述从中模拟数据的模型的曲线。

    Simulated predictor and response data for model 'Set 1'

+0

这正是我所需要的.U救了我! – Eddie 2012-02-22 18:26:11