使用VGLM伪R^2

Question

我被要求为我在R中运行的回归提供伪R2值。我使用vglm命令对4类患者生活质量数字运行回归，患者肺动脉应变，脉压和最大肺动脉压。

但是，vglm仅提供残余偏差并且不提供零偏差。我从以前的帖子中查到的伪R2，拟合度测量都需要无效偏差。有没有办法获得或计算与vglm的无效偏差？ https://www.stat.auckland.ac.nz/~yee/VGAM/doc/VGAMrefcard.pdf

谢谢， 肖纳

> summary(fit) 

Call: 
vglm(formula = resp_var ~ Strain + PP + Pmax, family = multinomial, data = DF) 

Pearson Residuals: 
        Min  1Q Median  3Q Max 
log(mu[,1]/mu[,4]) -3.8849 -0.45457 0.237302 0.456967 4.9936 
log(mu[,2]/mu[,4]) -2.7173 -0.45528 -0.255479 0.601062 2.1737 
log(mu[,3]/mu[,4]) -2.7175 -0.15174 -0.096334 -0.041271 6.6040 

Coefficients: 
       Estimate Std. Error z value 
(Intercept):1 -18.37926 1609.958 -0.0114160 
(Intercept):2 -22.19514 1609.958 -0.0137862 
(Intercept):3 -24.72028 1609.958 -0.0153546 
Strain:1  286.66750 3989.303 0.0718590 
Strain:2  285.95551 3989.303 0.0716806 
Strain:3  284.37742 3989.303 0.0712850 
PP:1   -0.87220  67.951 -0.0128356 
PP:2   -0.77133  67.951 -0.0113512 
PP:3   -0.93263  67.951 -0.0137250 
Pmax:1   0.15845  28.841 0.0054941 
Pmax:2   0.17548  28.841 0.0060843 
Pmax:3   0.28930  28.841 0.0100307 

Number of linear predictors: 3 

Names of linear predictors: 
log(mu[,1]/mu[,4]), log(mu[,2]/mu[,4]), log(mu[,3]/mu[,4]) 

Dispersion Parameter for multinomial family: 1 

Residual deviance: 102.8416 on 222 degrees of freedom 

Log-likelihood: -51.42081 on 222 degrees of freedom 

Number of iterations: 19 

Answer 1

空偏差是空（截距只）模型，例如残余偏差resp_var ~ 1。因此，要计算麦克法登的伪R平方：

# generate sample data 
set.seed(101) 
resp_var = t(rmultinom(100, 1, rep(0.25, 4))) 
DF = data.frame(Strain = rnorm(100), PP = rnorm(100), Pmax = rnorm(100)) 

# fit models 
m1 <- vglm(formula = resp_var ~ Strain + PP + Pmax, family = multinomial, data = DF) 
null <- vglm(formula = resp_var ~ 1, family = multinomial, data = DF) 

# McFadden's pseudo-R^2 
print(pseudo_R2 <- 1 - deviance(m1)/deviance(null)) 

# [1] 0.05300655 

此一致，例如，与mlogit：

library(mlogit) 
DF$cat_resp = factor(apply(resp_var, 1, function(.) which(. == 1))) 
mDf <- mlogit.data(DF, choice = "cat_resp", shape = "wide") 
summary(mlogit(cat_resp ~ 1|Strain + PP + Pmax, data = mDf)) 

# Call: 
# mlogit(formula = cat_resp ~ 1 | Strain + PP + Pmax, data = mDf, 
# method = "nr", print.level = 0) 
# ... 
# ... 
# Log-Likelihood: -130.03 
# McFadden R^2: 0.053007 
# Likelihood ratio test : chisq = 14.557 (p.value = 0.10386)