繪製適合R的GLM的輪廓偏差

我希望能夠繪製使用R中的函數glm()擬合的參數估計的輪廓偏差。輪廓偏差是對於參數估計的不同值的偏差函數在估計所有其他參數之後，有問題。我需要繪製擬合參數周圍幾個值的偏差，以檢查二次偏差函數的假設。繪製適合R的GLM的輪廓偏差

我的模型預測違法者的重新犯罪。該公式的格式爲： reconv ~ [other variables] + sex，其中reconv是二進制是/否因子，而sex是二進制男性/女性因子。我想繪製性別=女性（性別=男性是參考水平）估計的參數的配置偏差。

glm()函數估計參數爲-0.22，std誤差爲0.12。

[我在問這個問題，因爲沒有找到我能找到的答案，但是我解決了這個問題，並且希望發佈一個解決方案以供其他人使用。當然，歡迎提供其他幫助。 :-)]

來源

2012-03-20 MatW

查看Ioannis Kosmidis的profileModel package。他在R期刊上發表了一篇論文（將會出現R News）說明這個包裝：

Ioannis Kosmidis。 profilemodel R軟件包：爲具有線性預測器的模型分析目標。 R新聞，8（2）：12-18,2008年10月。

PDF是here（整個時事通訊）。

來源

2012-03-20 14:05:17

謝謝加文。這看起來就像我正在尋找的那種東西。我沒有意識到人們會很快回復。我的回答只是一些代碼，現在看起來有些多餘： -/ – MatW 2012-03-20 14:10:55

不，請張貼它---因爲許多不同方式的事情可以幫助未來毫無戒心的用戶。 – 2012-03-20 14:12:49

我剛剛發現我無法發佈七個小時，因爲我的聲譽不夠高。我會在稍後發佈。謝謝你的幫助。 – MatW 2012-03-20 14:15:11

見?profile.glm（和example("profile.glm")）在MASS包 - 我認爲它會做你想要的（這是默認不加載的一切，但它在?profile提到的，這可能是第一位你看起來......）（注意，剖面一般以符號平方根標度繪製，因此真正的二次剖面將顯示爲一條直線。）

來源

2012-03-20 13:58:56

謝謝本。我確實嘗試過，但並不完全理解情節告訴我的情況，因爲我期待的是一個倒二次形而不是一條直線。現在更有意義 - 謝謝。 – MatW 2012-03-20 14:05:52

好，夠公平的。將來，將這些細節添加到您的問題中會是明智的（例如，「我發現'profile.glm'，但它似乎沒有爲我的問題給出明智的答案」）。 – 2012-03-20 14:32:45

+1這對'profile.glm（）'@BenBolker;像上面這樣的事情發生在某個地方（簡歷？），而在回到日常工作之前，我花了幾分鐘時間用'profile.glm（）'。我錯過了「簽名平方根」位。 – 2012-03-20 15:01:31

我發現這樣做的方式涉及到使用offset（）函數（詳見Pawitan，Y。（2001）'In All Likelihood'p172）。由@BenBolker和@GavinSimpson給出的答案比這更好，因爲他們引用的軟件包將完成所有這些工作，還有更多。我發佈這個是因爲它以另一種方式圍繞它，還有，「手動」繪製東西有時對於學習來說很好。它教會了我很多。

sexi <- as.numeric(data.frame$sex)-1  #recode a factor as 0/1 numeric 

beta <- numeric(60)    #Set up vector to Store the betas 
deviance <- numeric(60)   #Set up vector to Store the deviances 

for (i in 1:60){ 

    beta[i] <- 0.5 - (0.01*i) 
    #A vector of values either side of the fitted MLE (in this case -0.22) 

    mod <- update(model, 
        .~. - sex    #Get rid of the fitted variable 
        + offset( I(sexi*beta[i]) ) #Replace with offset term. 
       ) 
    deviance[i] <- mod$deviance      #Store i'th deviance 
} 

best <- which.min(deviance)     
#Find the index of best deviance. Should be the fitted value from the model 

deviance0 <- deviance - deviance[best]   
#Scale deviance to zero by subtracting best deviance 

betahat <- beta[best] #Store best beta. Should be the fitted value. 
stderror <- 0.12187  #Store the std error of sex, found in summary(model) 

quadratic <- ((beta-betahat)^2)*(1/(stderror^2)) 
#Quadratic reference function to check quadratic assumption against 

x11()          
plot(beta,deviance0,type="l",xlab="Beta(sex)",ylim=c(0,4))  
lines(beta,quadratic,lty=2,col=3)   #Add quadratic reference line 
abline(3.84,0,lty=3)    #Add line at Deviance = 3.84

來源

2012-03-23 16:10:58 MatW

繪製適合R的GLM的輪廓偏差

回答

相關問題