Question sur nls fit dans R - pourquoi est-ce un ajustement si étrange?
J'essaie d'effectuer un ajustement non linéaire à certaines données simples (rendement du maïs par année). C'est assez simple de le faire avec lm dans R, mais certaines des données conviendraient mieux s'il y avait une courbe autorisée, quelque chose de l'ordre de l'année ^ 1,5 ou plus.
x <- c(1979L, 1980L, 1981L, 1982L, 1983L, 1984L, 1985L, 1986L, 1987L,
1988L, 1989L, 1990L, 1991L, 1992L, 1993L, 1994L, 1995L, 1996L,
1997L, 1998L, 1999L, 2000L, 2001L, 2002L, 2003L, 2004L, 2005L,
2006L, 2007L, 2008L, 2009L, 2010L, 2011L, 2012L, 2013L, 2015L,
2016L, 2017L, 2018L, 2019L)
y <- c(47.3, 25.4, 39, 56.4, 41.4, 56.1, 60.3, 58, 64, 35, 56, 54,
37, 80, 59, 88, 55, 87, 90, 99, 93, 90.4, 80.7, 35, 80.2, 104.9,
59.9, 43.5, 97.9, 106, 132, 121.7, 120.1, 63.9, 142.5, 129.9,
114.8, 122.1, 164.3, 133.9)
yield_model <- nls(y ~ x^a,start=list(a = 1))
plot(x,y)
lines(x,predict(yield_model),lty=2,col="red",lwd=3)
> yield_model2
Nonlinear regression model
model: y ~ x^a
data: parent.frame()
a
0.5778
residual sum-of-squares: 46984
Number of iterations to convergence: 8
Achieved convergence tolerance: 7.566e-09
Pourquoi les nls s'adaptent-ils si mal (visibles si vous les tracez)? Est-ce que j'ai fait quelque chose de mal? Vous pouvez imaginer qu'une légère courbe dans un ajustement aux données serait meilleure, avec une tendance. C'est comme si nls avait supprimé la tendance ou quelque chose du genre. Toute aide est la bienvenue.
Réponses
Deux options. Comme mentionné par @RuiBarradas, le problème est la spécification du modèle. Vous pouvez définir vos valeurs de départ en utilisant un lm()
comme ceci:
#Define initial values
mod <- lm(y~x)
#nls model
yield_model <- nls(y ~ a+x^b,
start=list(a = mod$coefficients[1],b=mod$coefficients[2]))
#Plot
plot(x,y)
lines(x,predict(yield_model),lty=2,col="red",lwd=3)
Production:

Ou essayez une autre approche comme loess
:
library(ggplot2)
#Data
df <- data.frame(x=x,y=y)
#Plot
ggplot(df,aes(x=x,y=y))+
geom_point()+
stat_smooth(se=F)
Production:

L'ajustement oublie le terme constant, l'ordonnée à l'origine. Contrairement aux autres fonctions de modélisation, nls
nécessite une interception explicite.
Ci-dessous, je monte également un modèle linéaire avec lm
, à titre de comparaison.
df1 <- data.frame(x, y)
yield_model <- nls(y ~ k + x^a, data = df1, start=list(k = 0, a = 1))
yield_model2 <- lm(y ~ x, df1)
summary(yield_model)
summary(yield_model2)
plot(x, y)
lines(x, predict(yield_model), lty = "dashed", col = "red", lwd = 3)
lines(x, predict(yield_model2), lty = "dotted", col = "blue", lwd = 3)

Comme vous pouvez le voir, les ajustements sont très proches les uns des autres. Mais ils ne sont pas égaux, à voir courir:
predict(yield_model) - predict(yield_model2)