Comment obtenir l'intervalle de confiance sur le changement du carré de la population

10

Pour un exemple simple, supposons qu'il existe deux modèles de régression linéaire

  • Modèle 1 a trois prédicteurs, x1a, x2betx2c
  • Le modèle 2 a trois prédicteurs du modèle 1 et deux prédicteurs supplémentaires x2aetx2b

Il existe une équation de régression de la population où la variance de la population expliquée est pour le modèle 1 et pour le modèle 2. La variance incrémentielle expliquée par le modèle 2 dans la population estρ(1)2ρ(2)2Δρ2=ρ(2)2ρ(1)2

Je souhaite obtenir des erreurs standard et des intervalles de confiance pour un estimateur de . Bien que l'exemple concerne respectivement 3 et 2 prédicteurs, mon intérêt de recherche concerne un large éventail de nombres différents de prédicteurs (par exemple, 5 et 30). Ma première pensée a été d'utiliser comme estimateur et de le bootstrap, mais je ne savais pas si cela Être approprié.Δρ2Δradj2=radj(2)2radj(1)2

Des questions

  • Est Δradj2 un estimateur raisonnable de Δρ2 ?
  • Comment obtenir un intervalle de confiance pour le changement du carré de la population (c.-à-d. Δρ2 )?
  • Le bootstrapping Δρ2 conviendrait-il pour le calcul de l'intervalle de confiance?

Toute référence à des simulations ou à la littérature publiée serait également la bienvenue.

Exemple de code

Si cela aide, j'ai créé un petit ensemble de données de simulation dans R qui pourrait être utilisé pour démontrer une réponse:

n <- 100
x <- data.frame(matrix(rnorm(n *5), ncol=5))
names(x) <- c('x1a', 'x1b', 'x1c', 'x2a', 'x2b')
beta <- c(1,2,3,1,2)
model2_rho_square <- .7
error_rho_square <- 1 - model2_rho_square
error_sd <- sqrt(error_rho_square / model2_rho_square* sum(beta^2))
model1_rho_square <- sum(beta[1:3]^2) / (sum(beta^2) + error_sd^2)
delta_rho_square <- model2_rho_square - model1_rho_square

x$y <- rnorm(n, beta[1] * x$x1a + beta[2] * x$x1b + beta[3] * x$x1c +
               beta[4] * x$x2a + beta[5] * x$x2b, error_sd)

c(delta_rho_square, model1_rho_square, model2_rho_square)
summary(lm(y~., data=x))$adj.r.square - 
        summary(lm(y~x1a + x1b + x1c, data=x))$adj.r.square

Raison de préoccupation avec bootstrap

J'ai exécuté un bootstrap sur certaines données avec environ 300 cas, et 5 prédicteurs dans le modèle simple et 30 prédicteurs dans le modèle complet. Bien que l'estimation de l'échantillon à l'aide de la différence r ajustée soit 0.116, l'intervalle de confiance boostrapped était pour la plupart plus grand IC95% (0,095 à 0,214) et la moyenne des bootstraps était loin de l'estimation de l'échantillon. La moyenne des échantillons boostés semble plutôt être centrée sur l'estimation de l'échantillon de la différence entre les carrés r dans l'échantillon. Ceci malgré le fait que j'utilisais les r-carrés ajustés de l'échantillon pour estimer la différence.

Fait intéressant, j'ai essayé une autre façon de calculer commeΔρ2

  1. calculer l'échantillon de changement de carré
  2. ajuster le changement de carré de l'échantillon en utilisant la formule standard de carré ajusté

Lorsqu'elle est appliquée aux données de l'échantillon, cela réduit l'estimation de à mais les intervalles de confiance semblent appropriés pour la méthode que j'ai mentionnée en premier, IC95% (.062, .179) avec une moyenne de .118.Δρ2.082

En gros, je crains que le bootstrap suppose que l'échantillon est la population, et donc les estimations que la réduction pour sur-ajustement peut ne pas fonctionner correctement.

Jeromy Anglim
la source
"Cependant, je crains que l'amorçage d'une telle valeur ajustée en fonction de la population puisse être problématique." -- Pourquoi?
Janvier
@Janvier J'ai édité la question et essayé d'exprimer ma préoccupation à propos du bootstrap avec un carré ajusté.
Jeromy Anglim
Quelle est la population R au carré ? J'ai jeté un œil à la définition donnée ici mais pour moi la variance n'a aucun sens car les ne sont pas distribués de façon identique. y iσy2yi
Stéphane Laurent
@ StéphaneLaurent c'est le pourcentage de variance expliqué dans la population par l'équation de régression de la population. Ou vous pouvez la définir asymptotiquement comme la proportion de variance expliquée dans votre échantillon lorsque votre taille d'échantillon approche de l'infini. Voir également cette réponse concernant les estimations non biaisées du carré de la population . Il est particulièrement pertinent en psychologie où nous sommes souvent plus intéressés par la vraie relation plutôt que d'appliquer réellement notre équation de prédiction estimée.
Jeromy Anglim
3
Un test F peut être considéré comme un test de l'hypothèse . Peut-il être utilisé pour dériver l'erreur standard et l'intervalle de confiance que vous recherchez? Δρ2=0
Maarten Buis

Réponses:

3

PopulationR2

Je suis tout d' abord essayer de comprendre la définition de la population R-carré .

Citant votre commentaire:

Ou vous pouvez la définir asymptotiquement comme la proportion de variance expliquée dans votre échantillon lorsque votre taille d'échantillon approche de l'infini.

Je pense que vous voulez dire que c'est la limite de l'échantillon lorsque l'on réplique le modèle infiniment de fois (avec les mêmes prédicteurs à chaque répétition). R2

Alors, quelle est la formule de la valeur asymptotique de l'échantillon ? Écrivez votre modèle linéaire comme dans https://stats.stackexchange.com/a/58133/8402 , et utilisez les mêmes notations que ce lien. On peut alors vérifier que l'échantillon va dans quand on réplique le modèle infiniment de fois.Y = μ + σ G R 2 p o p R 2 : = λR²Y=μ+σG
R2 Y=μ+σGpopR2:=λn+λY=μ+σG

Comme exemple:

> ## design of the simple regression model lm(y~x0)
> n0 <- 10
> sigma <- 1
> x0 <- rnorm(n0, 1:n0, sigma)
> a <- 1; b <- 2 # intercept and slope
> params <- c(a,b)
> X <- model.matrix(~x0)
> Mu <- (X%*%params)[,1]
> 
> ## replicate this experiment k times 
> k <- 200
> y <- rep(Mu,k) + rnorm(k*n0)
> # the R-squared is:
> summary(lm(y~rep(x0,k)))$r.squared 
[1] 0.971057
> 
> # theoretical asymptotic R-squared:
> lambda0 <- crossprod(Mu-mean(Mu))/sigma^2
> lambda0/(lambda0+n0)
          [,1]
[1,] 0.9722689
> 
> # other approximation of the asymptotic R-squared for simple linear regression:
> 1-sigma^2/var(y)
[1] 0.9721834

Population d'un sous-modèleR2

Supposons maintenant que le modèle est avec et considérons le sous-modèle . H1:μW1H0:μW0Y=μ+σGH1:μW1H0:μW0

Ensuite, j'ai dit ci-dessus que la population du modèle est où et et ensuite on a simplement .H 1 p o p R 2 1 : = λ 1R2H1 λ1= P Z 1 μ2popR12:=λ1n+λ1 Z1=[1]W1PZ1μ2=(μi-ˉμ)2λ1=PZ1μ2σ2Z1=[1]W1PZ1μ2=(μiμ¯)2

Définissez-vous maintenant la population du sous - modèle comme la valeur asymptotique du calculée par rapport au modèle mais sous l'hypothèse de distribution du modèle ? La valeur asymptotique (s'il y en a une) semble plus difficile à trouver.H 0 R 2 H 0 H 1R2 H0R2H0H1

Stéphane Laurent
la source
Merci Stéphane. Je vais devoir réfléchir à ce que vous dites. En ce qui concerne votre question. Je suppose que le véritable processus de génération de données n'est pas connu mais qu'il est le même pour les deux modèles, mais qu'il existe une véritable proportion de variance expliquée par la régression linéaire dans le modèle 1 et le modèle 2.
Jeromy Anglim
La formule @JeromyAnglim (A3) de cet article est un cas particulier de ma formule pour le modèle ANOVA unidirectionnel. Ma formule devrait donc être la définition générale de la population , mais ce n'est pas ce que vous utilisez dans votre PO. R2
Stéphane Laurent
1
@JeromyAnglim L'étude de cet article semble être proche de ce que vous recherchez (avec des prédicteurs aléatoires).
Stéphane Laurent
Merci. Le document d'Algina, Keselman et Penfield semble très utile. J'ai ajouté quelques commentaires à ma réponse à ce sujet.
Jeromy Anglim
@JeromyAnglim Alors, quelle est l'hypothèse concernant les prédicteurs? Ils sont générés selon une distribution gaussienne multivariée?
Stéphane Laurent
1

Plutôt que de répondre à la question que vous avez posée, je vais vous demander pourquoi vous posez cette question. Je suppose que vous voulez savoir si

mod.small <- lm(y ~ x1a + x1b + x1c, data=x)

est au moins aussi bon que

mod.large <- lm(y ~ ., data=x)

à expliquer y. Étant donné que ces modèles sont imbriqués, la façon évidente de répondre à cette question semble être d'exécuter une analyse de variance en les comparant, de la même manière que vous pourriez exécuter une analyse de la déviance pour deux GLM, comme

anova(mod.small, mod.large)

Ensuite, vous pouvez utiliser l'échantillon d'amélioration du carré R entre les modèles comme meilleure estimation de ce que serait l'amélioration de l'ajustement dans la population, en supposant toujours que vous pouvez donner un sens à la population R au carré. Personnellement, je ne suis pas sûr de pouvoir le faire, mais cela n'a pas d'importance dans les deux cas.

Plus généralement, si vous êtes intéressé par les quantités de population, vous êtes probablement intéressé par la généralisation, donc une mesure d'ajustement d'échantillon n'est pas tout à fait ce que vous voulez, même si elle est «corrigée». Par exemple, la validation croisée d'une certaine quantité qui estime le type et la quantité d'erreurs réelles que vous pourriez vous attendre à faire à partir d'un échantillon, comme MSE, semblerait obtenir ce que vous voulez.

Mais il est fort possible que je manque quelque chose ici ...

conjugateprior
la source
J'apprécie votre réponse, et ce pourrait bien être un bon conseil pour les autres. Mais mon contexte de recherche signifie que je suis légitimement intéressé par la place delta-rho. Alors que la plupart des statisticiens sont souvent plus préoccupés par l'utilité prédictive d'un modèle (par exemple, le carré r delta validé croisé), je suis un psychologue et je suis spécifiquement intéressé par la propriété de la population. De plus, je ne suis pas intéressé par la signification statistique de l'amélioration. Je m'intéresse à la taille de l'amélioration. Et je trouve que delta-r-square est une métrique utile pour indexer cette taille d'amélioration.
Jeromy Anglim
En ce qui concerne le MSE, différentes études en psychologie utilisent des mesures sur des paramètres très différents. Ainsi, il existe une attirance, à tort ou à raison, pour les mesures normalisées telles que le r-carré.
Jeromy Anglim
Assez juste, en particulier sur MSE. Je reste un peu confus par l'intérêt pour le bootstrap et l'inférence de la population mais le manque d'intérêt pour les tests car, peut-être naïvement, ces préoccupations semblent être des préoccupations équivalentes traitées différemment. J'ai également de la difficulté à distinguer étroitement la prédiction de l'échantillon de l'inférence à une population, mais c'est probablement le bayésianisme avant le café (où la prédiction n'est qu'un autre problème d'inférence de la population) se met en travers.
conjugateprior
J'ai peut-être parlé un peu rapidement. Dans mon contexte de recherche, il y a souvent beaucoup de preuves que le delta-rho-carré est supérieur à zéro. La question d'intérêt est quel est le degré d'augmentation. C'est-à-dire qu'il s'agit d'une augmentation triviale ou d'une augmentation théoriquement significative. Ainsi, la confiance ou des intervalles crédibles me donnent une estimation de l'incertitude autour de cette augmentation. Je n'ai pas encore concilié ce que je fais ici avec ma compréhension des statistiques bayésiennes, mais j'aimerais bien.
Jeromy Anglim
1

Les éléments suivants représentent quelques possibilités de calcul des intervalles de confiance sur .ρ2

Bootstrap carré double ajusté

Ma meilleure supposition actuelle sur une réponse est de faire un bootstrap r-square à double ajustement. J'ai implémenté la technique. Elle implique les éléments suivants:

  • Générez un ensemble d'échantillons d'amorçage à partir des données actuelles.
  • Pour chaque échantillon amorcé:
    • calculer le premier carré r ajusté pour les deux modèles
    • calculer le deuxième carré r ajusté sur les valeurs de carré r ajustées de l'étape précédente
    • Soustrayez model2 des secondes valeurs de r ajustées du modèle1 pour obtenir une estimation de .Δρ2

Le raisonnement est que le premier carré r ajusté supprime le biais introduit par le bootrapping (c.-à-d. Que le bootstrapping suppose que le carré r de l'échantillon est le carré r de la population). Le deuxième carré r ajusté effectue la correction standard qui est appliquée à un échantillon normal pour estimer le carré r de la population.

À ce stade, tout ce que je peux voir, c'est que l'application de cet algorithme génère des estimations qui semblent à peu près correctes (c'est-à-dire que le theta_hat moyen dans le bootstrap est très proche de l'exemple theta_hat). L'erreur standard correspond à mon intuition. Je n'ai pas encore testé s'il fournit une couverture fréquentiste appropriée là où le processus de génération de données est connu, et je ne suis pas non plus entièrement sûr à ce stade comment l'argument pourrait être justifié à partir des premiers principes

Si quelqu'un voit des raisons pour lesquelles cette approche serait problématique, je serais reconnaissant d'en entendre parler.

Simulation par Algina et al

Stéphane a mentionné l'article d'Algina, Keselman et Penfield. Ils ont effectué une étude de simulation pour examiner la couverture de l'intervalle de confiance à 95% des méthodes d'amorçage et asymptotiques pour estimer . Leurs méthodes d'amorçage n'impliquaient qu'une seule application du carré r ajusté, plutôt que le double ajustement du carré r que je mentionne ci-dessus. Ils ont constaté que les estimations bootstrap ne fournissaient une bonne couverture que lorsque le nombre de prédicteurs supplémentaires dans le modèle complet était de un ou peut-être deux. C'est mon hypothèse que c'est parce que plus le nombre de prédicteurs augmente, plus la différence entre le bootstrap r et carré ajusté simple et double augmente.Δρ2

Smithson (2001) sur l'utilisation du paramètre de non-centralité

Smithson (2001) discute du calcul des intervalles de confiance pour le partiel en fonction du paramètre de non-centralité. Voir notamment pages 615 et 616. Il suggère qu '"il est simple de construire un IC pour et partiel mais pas pour la corrélation semi-partisane au carré". (p.615)f 2 R 2R2f2R2

Références

  • Algina, J., Keselman, HJ et Penfield, RD Intervalles de confiance pour le coefficient de corrélation semipartiale multiple au carré. PDF
  • Smithson, M. (2001). Intervalles de confiance corrects pour diverses tailles et paramètres d'effet de régression: l'importance des distributions non centrales dans le calcul des intervalles. Mesures éducatives et psychologiques, 61 (4), 605-632.
Jeromy Anglim
la source
1
Il semble que personne ici (y compris vous) ne connaisse la définition de votre population au carré R. Par conséquent, à mon humble avis, il s'agit d'une approche très problématique.
Stéphane Laurent
@ StéphaneLaurent Merci pour cela. J'avoue que jusqu'à présent, je n'ai pas considéré le carré de population comme une propriété de discorde. Par exemple, je pourrais proposer un processus de génération de données et il y aurait un carré r qui est approché lorsque ma taille d'échantillon de simulation approche de l'infini. Et de même, je suppose qu'il existe un processus de génération de données pour mes données, et donc s'il était possible d'obtenir un échantillon infini, je pourrais calculer le vrai carré de la population.
Jeromy Anglim
Oui, mais j'ai l'impression que vous supposez également un processus de génération pour les prédicteurs. Je ne peux pas comprendre comment cela pourrait avoir un sens pour un modèle linéaire général.
Stéphane Laurent