R语言中用于计算Rsquare的包rsq

路思源

2023-12-01

理论介绍

线性模型情形

R-squared（值范围 0-1）描述的输入变量对输出变量的解释程度。在单变量线性回归中R-squared 越大，说明拟合程度越好。

数学表达式：

$R^2=SSR/TSS=1-RSS/TSS$
其中：

TSS 是执行回归分析前，响应变量固有的方差
RSS 是残差平方和（就是回归模型不能解释的方差）
SSR 回归模型可以解释的方差

然而，只要增加了更多的变量，无论增加的变量是否和输出变量存在关系，则R-squared 要么保持不变，要么增加。

所以需要 adjusted R-squared （范围 $(-\infty,1]$ ）,它会对那些增加的且不会改善模型效果的变量增加一个惩罚项。

数学表达式：

$R^2(adj)=1-(1-R^2)\times\frac{n-1}{n-p-1}$

(p 为变量个数，n 为样本个数)

另外一个表达是：

$R^2(adj)=1-\frac{RSS/(n-p-1)}{TSS/(n-1)}$

结论，如果单变量线性回归，则使用 R-squared评估，多变量，则使用adjusted R-squared。

在单变量线性回归中，R-squared和adjusted R-squared是一致的。

另外，如果增加更多无意义的变量，则 R-squared 和adjusted R-squared之间的差距会越来越大，Adjusted R-squared会下降。但是如果加入的特征值是显著的，则adjusted R-squared 也会上升。
参考来源：https://zhuanlan.zhihu.com/p/340105912

R语言中的req包可以对线性模型和广义线性模型的多种形式的 $R^2$ 进行计算

广义线性模型情形

针对广义线性模型R-squared通常定义为解释的不确定性的比例（the proportionate reduction in uncertainty, measured by KL divergence, due to the inclusion of regressors.）
一种基于KL散度的R-squared定义是：
$R_{KL}^2=1-\frac{K(y,\hat{\mu})}{K(y,\hat{\mu}_0)}$ ,
其中，
$K(y,\hat{\mu})=\sum_{i=1}^n[\log f_y(y_i)-\log f_{\hat{\mu}}(y_i)]$ ,
$\hat{\mu}_i=\exp(x_i \hat{\beta})$ 表示模型估计的均值
$\hat{\mu}_0$ 是对 $\mu$ 直接估计的极大似然估计，往往是 $y$ 的均值

参考文献：
A. C. Cameron, F. A. G. Windmeijer. An R-squared measure of goodness of fit for some common nonlinear regression models. Journal of Econometrics. 1997, 77: 329-342.

函数介绍

下面的函数介绍来源于：https://cran.r-project.org/web/packages/rsq/rsq.pdf

rsq()的介绍

rsq(fitObj,adj=FALSE,type=c('v','kl','sse','lr','n'))

参数介绍：

fitObj：一个“lm”,“glm”, “merMod”, “lmerMod” 或者类 “lme”; 通常是 lm, glm, glm.nb, lmer , glmer , lme4中的glmer ，nlme中 lme的一个结果。
adj：这是一个逻辑值，如果是TRUE，则计算adjusted R^2
type: 表示R-squared的类型 (仅对广义线性模型有用):
’v’ (default) – variance-function-based (Zhang, 2016), calling rsq.v;
’kl’ – KL-divergence-based (Cameron and Windmeijer, 1997), calling rsq.kl;
’sse’ – SSE-based (Efron, 1978), calling rsq.sse;
’lr’ – likelihood-ratio-based (Maddala, 1983; Cox and Snell, 1989; Magee,
1990), calling rsq.lr;
’n’ – corrected version of ’lr’ (Nagelkerke, 1991), calling rsq.n

输出参数：
该函数除了返回R-square的值，针对（广义）线性混合模型还有下列输出：

R_M^2: 模型总共可以解释的变异性的比例，包含固定效应以及随机效应.
R_F^2：可以被固定因子解释的变异性所占的比例
R_R^2：可以被随机效应因子解释的变异性所占的比例

实例

library(rsq)
data(hcrabs)
attach(hcrabs)
y <- ifelse(num.satellites>0,1,0)
bnfit <- glm(y~color+spine+width+weight,family=binomial)
rsq(bnfit)
# [1] 0.2171238
rsq(bnfit,adj=TRUE)
# [1] 0.1839109
quasibn <- glm(y~color+spine+width+weight,family=quasibinomial)
rsq(quasibn)
# [1] 0.2171238
rsq(quasibn,adj=TRUE)
# [1] 0.1839109
psfit <- glm(num.satellites~color+spine+width+weight,family=poisson)
rsq(psfit)
# [1] 0.1172267
rsq(psfit,adj=TRUE)
# [1] 0.07977572
quasips <- glm(num.satellites~color+spine+width+weight,family=quasipoisson)
rsq(quasips)
# [1] 0.1172267
rsq(quasips,adj=TRUE)
# [1] 0.07977572



# Linear mixed models
require(lme4)
lmm1 <- lmer(Reaction~Days+(Days|Subject),data=sleepstudy)
rsq(lmm1)
# $model
# [1] 0.8003832
# 
# $fixed
# [1] 0.2864714
# 
# $random
# [1] 0.5139119
rsq.lmm(lmm1)
# $model
# [1] 0.8003832
# 
# $fixed
# [1] 0.2864714
# 
# $random
# [1] 0.5139119

rsq.partial()的介绍

该函数用于计算aka partial $R^2$ , 可用于线性模型和广义线性模型

rsq.partial(objF,objR=NULL,adj=FALSE,type=c('v','kl','sse','lr','n'))

参数介绍：

objF：是一个 “lm” 或 "glm"对象, 也是lm, glm, 或者glm.nb 的用于拟合全模型的结果
objR：是一个 “lm” 或 "glm"对象, 也是lm, glm, 或者glm.nb 的用于拟合缩减模型的结果
adj: logical; if TRUE, calculate the adjusted partial R^2.
type: R-squared的类型:
’v’ (default) – variance-function-based (Zhang, 2016), calling rsq.v;
’kl’ – KL-divergence-based (Cameron and Windmeijer, 1997), calling rsq.kl;
’sse’ – SSE-based (Efron, 1978), calling rsq.sse;
’lr’ – likelihood-ratio-based (Maddala, 1983; Cox and Snell, 1989; Magee,1990), calling rsq.lr;
’n’ – corrected version of ’lr’ (Nagelkerke, 1991), calling rsq.n.

输出参数：
返回值包含adjustment 和partial.rsq. 当objR 为 NULL时, variable.full和 variable.reduced将被返回; 否则返回的是variable。

adjustment：逻辑值; 如果为TRUE，则计算adjusted partial R^2.
variable.full ：全模型中全部的协变量
variable.reduced：缩减模型中的全部协变量.
variable：全模型中全部的协变量.
partial.rsq： partial $R^2$ 或者 adjusted partial $R^2$ .

pcor()函数介绍

该函数用于计算线性模型和广义线性模型的偏相关系数

pcor(objF,objR=NULL,adj=FALSE,type=c('v','kl','sse','lr','n'))

参数介绍：

objF：是一个 “lm” 或 "glm"对象, 也是lm, glm, 或者glm.nb 的用于拟合全模型的结果
objR：是一个 “lm” 或 "glm"对象, 也是lm, glm, 或者glm.nb 的用于拟合缩减模型的结果
adj: logical; if TRUE, calculate the adjusted partial R^2.
type: R-squared的类型:
’v’ (default) – variance-function-based (Zhang, 2016), calling rsq.v;
’kl’ – KL-divergence-based (Cameron and Windmeijer, 1997), calling rsq.kl;
’sse’ – SSE-based (Efron, 1978), calling rsq.sse;
’lr’ – likelihood-ratio-based (Maddala, 1983; Cox and Snell, 1989; Magee,1990), calling rsq.lr;
’n’ – corrected version of ’lr’ (Nagelkerke, 1991), calling rsq.n.

注意：
当缩减模型的拟合对象缺失时，对于每个变量的偏相关系数都会被计算（除了多于两个水平的因子）。

vresidual()函数介绍

用于计算variance-function-based的残差, 这是用来计算variance-function-based R-squared.

vresidual(y,yfit,family=binomial(),variance=NULL)

参数介绍：

y：观测值的向量
yfit：拟合值的向量
family：分布族
variance：方差函数 (specified by family by default).

注意：
残差的计算依赖于方差函数, 对于quasi 模型需要很好的定义. 当方差函数是常数或者线性函数时，这会导出经典的残差。注意只有方差函数需要被设定，通过 “family”" 或者 “variance”。

实例

data(hcrabs)
attach(hcrabs)
y <- ifelse(num.satellites>0,1,0)
bnfit <- glm(y~color+spine+width+weight,family="binomial")
vresidual(y,bnfit$fitted.values,family="binomial")
# Effectiveness of Bycycle Safety Helmets in Thompson et al. (1989)
y <- matrix(c(17,218,233,758),2,2)
x <- factor(c("yes","no"))
tbn <- glm(y~x,family="binomial")
yfit <- cbind(tbn$fitted.values, 1-tbn$fitted.values)
vr0 <- vresidual(matrix(0,2,1),yfit[,1],family="binomial")
vr1 <- vresidual(matrix(1,2,1),yfit[,2],family="binomial")
y[,1]*vr0+y[,2]*vr1

R语言中用于计算Rsquare的包rsq

理论介绍

线性模型情形

广义线性模型情形

函数介绍

rsq()的介绍

实例

rsq.partial()的介绍

相关文献

pcor()函数介绍

vresidual()函数介绍

实例

相关阅读

相关文章

相关问答

相关文档