R将nnls转换为lm

Question

我正在寻找一种在积极约束下进行线性回归的方法,因此遇到了nnls方法.但是我想知道如何从nnls获得与lm对象提供的统计信息相同的统计信息.更具体地说,R平方,akaike标准,p值.

library(arm)
library(nnls)


data = runif(100*4, min = -1, max = 1)
data = matrix(data, ncol = 4)
colnames(data) = c("y", "x1", "x2", "x3")
data = as.data.frame(data)
data$x1 = -data$y

A = as.matrix(data[,c("x1", "x2", "x3")])
b = data$y

test = nnls(A,b)
print(test)

有没有办法在lm框架中重新估计,使用偏移和修复系数不起作用...有没有办法获得这些统计数据？或者在系数的正约束下创建lm对象的另一种方法？

谢谢罗曼.

Answer 1

你提议做的是一个非常糟糕的想法,以至于我不愿意告诉你如何去做.原因是对于OLS,假设残差正态分布且方差不变,则参数估计遵循多变量t分布,我们可以通常的方式计算置信限和p值.

但是,如果我们对相同的数据执行NNLS,则残差通常不会被分配,并且用于计算p值等的标准技术将产生垃圾.有一些方法可以估算NNLS拟合参数的置信限(例如参见此参考),但它们是近似的,通常依赖于对数据集的相当严格的假设.

在另一方面,这将是很好,如果一些为更基本的功能lm对象,如predict(...),coeff(...),residuals(...),等还曾为一个NNLS配合的结果.因此,实现这一点的一种方法是使用nls(...):仅仅因为模型在参数中是线性的并不意味着您不能使用非线性最小二乘来查找参数.nls(...)如果使用port算法,则可以选择设置参数的下限(和上限).

set.seed(1)   # for reproducible example
data   <- as.data.frame(matrix(runif(1e4, min = -1, max = 1),nc=4))
colnames(data) <-c("y", "x1", "x2", "x3")
data$y <- with(data,-10*x1+x2 + rnorm(2500)) 

A <- as.matrix(data[,c("x1", "x2", "x3")])
b <- data$y
test <- nnls(A,b)
test
# Nonnegative least squares model
# x estimates: 0 1.142601 0 
# residual sum-of-squares: 88391
# reason terminated: The solution has been computed sucessfully.

fit <- nls(y~b.1*x1+b.2*x2+b.3*x3,data,algorithm="port",lower=c(0,0,0))
fit
# Nonlinear regression model
#   model: y ~ b.1 * x1 + b.2 * x2 + b.3 * x3
#    data: data
#   b.1   b.2   b.3 
# 0.000 1.143 0.000 
#  residual sum-of-squares: 88391

正如你所看到的,使用的结果nnls(...),并使用的结果,nls(...)与lower-c(0,0,0)是相同的.但是nls(...)生成一个nls对象,它支持(大多数)与lm对象相同的方法.所以,你可以写precict(fit),coef(fit),residuals(fit),AIC(fit)等你也可以写summary(fit)和confint(fit) 但要注意:你得到的价值是没有意义的!

为了说明关于残差的点,我们将OLS拟合的残差与该数据进行比较,并将NNLS拟合的残差进行比较.

par(mfrow=c(1,2),mar=c(3,4,1,1))
qqnorm(residuals(lm(y~.,data)),main="OLS"); qqline(residuals(lm(y~.,data)))
qqnorm(residuals(fit),main="NNLS"); qqline(residuals(fit)) 

在该数据集中,变异性的随机部分y是设计的N(0,1),因此OLS拟合的残差(左边的QQ图)是正常的.但使用NNLS拟合的同一数据集的残差并不是很正常.这是因为yon 的真正依赖性x1是-10,但NNLS拟合迫使它为0.因此,非常大的残差(正面和负面)的比例远高于正态分布预期的比例.