比较完整模型与nls的空模型。

Question

我想比较一个完整的非线性模型和一个空模型。例如，考虑到通常使用glm的方法，这可能吗？在我的情况下。

#Packages
library(minpack.lm)

# Data set - Diameter in function of Feature and Age
Feature<-sort(rep(c("A","B"),22))
Age<-c(60,72,88,96,27,
36,48,60,72,88,96,27,36,48,60,72,
88,96,27,36,48,60,27,27,36,48,60,
72,88,96,27,36,48,60,72,88,96,27,
36,48,60,72,88,96)
Diameter<-c(13.9,16.2,
19.1,19.3,4.7,6.7,9.6,11.2,13.1,15.3,
15.4,5.4,7,9.9,11.7,13.4,16.1,16.2,
5.9,8.3,12.3,14.5,2.3,5.2,6.2,8.6,9.3,
11.3,15.1,15.5,5,7,7.9,8.4,10.5,14,14,
4.1,4.9,6,6.7,7.7,8,8.2)
d<-dados <- data.frame(Feature,Age,Diameter)
str(d)

# Complet model
e1<- Diameter ~ a1 * Age^a2 
#Algoritm Levenberg-Marquardt
m1 <-  nlsLM(e1, data = d,
     start = list(a1 = 0.1, a2 = 10),
     control = nls.control(maxiter = 1000))

#Null model

e2<- Diameter ~ 1 
#Algoritm Levenberg-Marquardt
m0 <-  nlsLM(e1, data = d,
     control = nls.control(maxiter = 1000))
Warning message:
In nlsLM(e1, data = d, control = nls.control(maxiter = 1000)) :
  No starting values specified for some parameters.
Initializing ‘a1’, ‘a2’ to '1.'.
Consider specifying 'start' or using a selfStart model

没有工作，我的最终目标是：

anova(m1,m0)

这在非线性宇宙中可能吗？

Answer 1

是的，有一个 anova.nls 定义于R

methods("anova")
## [1] anova.glm*      anova.glmlist*  anova.lm*       anova.lmlist*  
## [5] anova.loess*    anova.mlm*      anova.nls*      anova.quantmod*
## see '?methods' for accessing help and source code

但是，问题中的代码有一些问题。

• m1应该使用e2。按照写法，m0和m1都是同一个模型。
• e2的指定是错误的。 在问题中，它被指定为好像是一个lm模型，而它应该被指定为一个ns模型。
• 我会小心使用nsLM。 人们经常使用它，因为它似乎收敛了，而ns不收敛，但这可能发生，因为它过早地收敛了，给出了一个偏离的答案。

试试这个。

o <- order(d$Age)

# Complet model
e1<- Diameter ~ Age^a2 
m1a <-  nls(e1, data = d[o, ], start = list(a2 = 10), alg = "plinear")

#Null model

ones <- rep(1, nrow(d))
e2 <- Diameter ~ a * ones
m0a <-  nls(e2, data = d[o, ], start = list(a = 1))

anova(m1a, m0a)
## Analysis of Variance Table
##
## Model 1: Diameter ~ Age^a2
## Model 2: Diameter ~ a * ones
##   Res.Df Res.Sum Sq Df  Sum Sq F value   Pr(>F)    
## 1     42     270.96                                
## 2     43     823.73 -1 -552.77   85.68 1.07e-11 ***
## ---
## Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

plot(Diameter ~ Age, d)
lines(fitted(m0a) ~ Age, d[o, ], col = "blue")
lines(fitted(m1a) ~ Age, d[o, ], col = "red")