Comment puis-je trouver une valeur de p de régression lisse spline / loess?

10

J'ai quelques variables et je suis intéressé à trouver des relations non linéaires entre elles. J'ai donc décidé d'adapter des splines ou des loess, et d'imprimer de jolis tracés (voir le code ci-dessous). Mais, je veux aussi avoir des statistiques qui me donnent une idée de la probabilité que la relation soit une question de hasard ... c'est-à-dire, j'ai besoin d'une valeur de p globale, comme je l'ai pour la régression linéaire par exemple. En d'autres termes, j'ai besoin de savoir si la courbe ajustée a un sens, car mon code adaptera une courbe à toutes les données.

x <- rnorm(1000)
y <- sin(x) + rnorm(1000, 0, 0.5)

cor.test(x,y)
plot(x, y, xlab = xlab, ylab = ylab)
spl1 <- smooth.spline(x, y, tol = 1e-6, df = 8)
lines(spl1, col = "green", lwd = 2)

spl2 <- loess(y ~ x)
x.pr <- seq(min(x), max(x), length.out = 100)
lines(x.pr, predict(spl2, x.pr), col = "blue", lwd = 2)
Curieuse
la source

Réponses:

8

La bibliothèque de splines a des fonctions bset nscela créera une base de splines à utiliser avec la lmfonction, alors vous pouvez adapter un modèle linéaire et un modèle comprenant des splines et utiliser la anovafonction pour effectuer le test de modèle complet et réduit pour voir si le modèle de spline s'adapte nettement mieux que le modèle linéaire.

Voici un exemple de code:

x <- rnorm(1000)
y <- sin(x) + rnorm(1000, 0, 0.5)

library(splines)

fit1 <- lm(y~x)
fit0 <- lm(y~1)
fit2 <- lm(y~bs(x,5))

anova(fit1,fit2)
anova(fit0,fit2)

plot(x,y, pch='.')
abline(fit1, col='red')
xx <- seq(min(x),max(x), length.out=250)
yy <- predict(fit2, data.frame(x=xx))
lines(xx,yy, col='blue')

Vous pouvez également utiliser la polyfonction pour effectuer un ajustement polynomial et tester les termes non linéaires comme test de courbure.

Pour l'ajustement loess c'est un peu plus compliqué. Il existe des estimations de degrés de liberté équivalents pour le paramètre de lissage du loess qui pourraient être utilisées avec les valeurs pour les modèles linéaire et loess à construire et à tester F. Je pense que les méthodes basées sur les tests d'amorçage et de permutation peuvent être plus intuitives.R2

Il existe des techniques pour calculer et tracer un intervalle de confiance pour un ajustement de loess (je pense qu'il peut y avoir un moyen intégré dans le package ggplot2), vous pouvez tracer la bande de confiance et voir si une ligne droite rentrerait dans la bande (ceci n'est pas une valeur de p, mais donne toujours un oui / non.

Vous pouvez ajuster un modèle linéaire et prendre les résidus et ajuster un modèle de lœss aux résidus comme réponse (et la variable d'intérêt comme prédicteur), si le vrai modèle est linéaire, cet ajustement doit être proche d'une ligne plate et réorganiser les points par rapport au prédicteur ne devrait pas faire de différence. Vous pouvez l'utiliser pour créer un test de permutation. Ajustez le loess, trouvez la valeur prédite la plus éloignée de 0, permutez maintenant aléatoirement les points et ajustez un nouveau loess et trouvez le point le plus éloigné de 0, répétez plusieurs fois, la valeur p est la proportion de valeurs permutées qui sont plus loin de 0 à la valeur d'origine.

Vous pouvez également envisager la validation croisée comme méthode de choix de la bande passante loess. Cela ne donne pas de valeur p, mais une bande passante infinie correspond à un modèle linéaire parfait, si la validation croisée suggère une très large bande passante, alors cela suggère qu'un modèle linéaire peut être raisonnable, si les bandes passantes plus élevées sont clairement inférieures à certaines les bandes passantes plus petites suggèrent alors une courbure définie et linéaire n'est pas suffisant.

Greg Snow
la source
Merci Greg! Je pense que le 1er paragraphe semble être la voie à suivre, sauf que je ne suis pas intéressé par la comparaison avec le modèle linéaire, juste pour voir si la spline l'explique ou non. Pourriez-vous s'il vous plaît fournir du code ou des pointeurs plus concrets sur la façon de tester la spline avec anova? J'ai regardé les fonctions bs et ns mais je ne suis pas si bon en statistique pour pouvoir l'inventer moi-même.
Curieux
Et oui, j'ai trouvé comment calculer pour loess (voir fibosworld.wordpress.com/2012/11/04/loess-regression-with-r ) mais je ne sais pas comment convertir en p- valeur ...R2R2
Curieux
J'ai ajouté un exemple de code pour l' anovaapproche avec des splines. Pour le test F de considérons que est le SSR divisé par SST et est SSE divisé par SST, donc le rapport est juste le SSR divisé par SSE (les 2 cas de SST s'annulent). Incluez les degrés de liberté et une petite algèbre et vous avez la statistique F pour la régression globale. R2R21-R2R21-R2
Greg Snow
Greg, merci! 1) Pourriez-vous expliquer ce qui se lm(y~bs(x,5))passe et pourquoi ce n'est pas le cas lm(y~I(bs(x,5)))? Je suis assez confus par cet appel car le résultat de bs (x, 5) n'est pas une variable ... 2) Est-ce que je comprends bien que la valeur de p que je recherche est le résultat de anova(fit0,fit2)?
Curieux
1
Il s'agit simplement d'une régression linéaire. Tout comme vous pouvez faire une régression linéaire avec , et pour ajuster une courbe basée sur un polynôme, la fonction crée des transformations de la variable et les transmet à qui fait la régression linéaire. XX2X3bsXlm
Greg Snow