Les statisticiens supposent-ils qu'on ne peut pas trop arroser une plante, ou est-ce que j'utilise simplement les mauvais termes de recherche pour la régression curviligne?

18

Presque tout ce que j'ai lu sur la régression linéaire et GLM se résume à ceci: où est une fonction non croissante ou non décroissante de et est le paramètre que vous estimer et tester des hypothèses sur. Il existe des dizaines de fonctions de liaison et de transformations de et pour faire de une fonction linéaire de .y=f(x,β)f(x,β)xβyxyf(x,β)

Maintenant, si vous supprimez l'exigence non croissante / non décroissante pour , je ne connais que deux choix pour ajuster un modèle paramétrique linéarisé: les fonctions trigonométriques et les polynômes. Les deux créent une dépendance artificielle entre chaque prédit et l'ensemble entier de , ce qui en fait un ajustement très non robuste, sauf s'il existe des raisons préalables de croire que vos données sont réellement générées par un processus cyclique ou polynomial.f(x,β)XyX

Ce n'est pas une sorte de cas de bord ésotérique. C'est la relation réelle et sensée entre l'eau et les rendements des cultures (une fois que les parcelles sont suffisamment profondes sous l'eau, les rendements des cultures commenceront à diminuer), ou entre les calories consommées au petit-déjeuner et les performances d'un quiz mathématique, ou le nombre de travailleurs dans une usine et le nombre de widgets qu'ils produisent ... en bref, presque tous les cas de la vie réelle pour lesquels des modèles linéaires sont utilisés, mais avec des données couvrant une gamme suffisamment large pour que vous puissiez dépasser les rendements décroissants en rendements négatifs.

J'ai essayé de chercher les termes «concave», «convexe», «curviligne», «non monotone», «baignoire», et j'oublie combien d'autres. Peu de questions pertinentes et encore moins de réponses utilisables. Donc, en termes pratiques, si vous aviez les données suivantes (code R, y est une fonction de la variable continue x et du groupe de variables discrètes):

updown<-data.frame(y=c(46.98,38.39,44.21,46.28,41.67,41.8,44.8,45.22,43.89,45.71,46.09,45.46,40.54,44.94,42.3,43.01,45.17,44.94,36.27,43.07,41.85,40.5,41.14,43.45,33.52,30.39,27.92,19.67,43.64,43.39,42.07,41.66,43.25,42.79,44.11,40.27,40.35,44.34,40.31,49.88,46.49,43.93,50.87,45.2,43.04,42.18,44.97,44.69,44.58,33.72,44.76,41.55,34.46,32.89,20.24,22,17.34,20.14,20.36,24.39,22.05,24.21,26.11,28.48,29.09,31.98,32.97,31.32,40.44,33.82,34.46,42.7,43.03,41.07,41.02,42.85,44.5,44.15,52.58,47.72,44.1,21.49,19.39,26.59,29.38,25.64,28.06,29.23,31.15,34.81,34.25,36,42.91,38.58,42.65,45.33,47.34,50.48,49.2,55.67,54.65,58.04,59.54,65.81,61.43,67.48,69.5,69.72,67.95,67.25,66.56,70.69,70.15,71.08,67.6,71.07,72.73,72.73,81.24,73.37,72.67,74.96,76.34,73.65,76.44,72.09,67.62,70.24,69.85,63.68,64.14,52.91,57.11,48.54,56.29,47.54,19.53,20.92,22.76,29.34,21.34,26.77,29.72,34.36,34.8,33.63,37.56,42.01,40.77,44.74,40.72,46.43,46.26,46.42,51.55,49.78,52.12,60.3,58.17,57,65.81,72.92,72.94,71.56,66.63,68.3,72.44,75.09,73.97,68.34,73.07,74.25,74.12,75.6,73.66,72.63,73.86,76.26,74.59,74.42,74.2,65,64.72,66.98,64.27,59.77,56.36,57.24,48.72,53.09,46.53),
                   x=c(216.37,226.13,237.03,255.17,270.86,287.45,300.52,314.44,325.61,341.12,354.88,365.68,379.77,393.5,410.02,420.88,436.31,450.84,466.95,477,491.89,509.27,521.86,531.53,548.11,563.43,575.43,590.34,213.33,228.99,240.07,250.4,269.75,283.33,294.67,310.44,325.36,340.48,355.66,370.43,377.58,394.32,413.22,428.23,436.41,455.58,465.63,475.51,493.44,505.4,521.42,536.82,550.57,563.17,575.2,592.27,86.15,91.09,97.83,103.39,107.37,114.78,119.9,124.39,131.63,134.49,142.83,147.26,152.2,160.9,163.75,172.29,173.62,179.3,184.82,191.46,197.53,201.89,204.71,214.12,215.06,88.34,109.18,122.12,133.19,148.02,158.72,172.93,189.23,204.04,219.36,229.58,247.49,258.23,273.3,292.69,300.47,314.36,325.65,345.21,356.19,367.29,389.87,397.74,411.46,423.04,444.23,452.41,465.43,484.51,497.33,507.98,522.96,537.37,553.79,566.08,581.91,595.84,610.7,624.04,637.53,649.98,663.43,681.67,698.1,709.79,718.33,734.81,751.93,761.37,775.12,790.15,803.39,818.64,833.71,847.81,88.09,105.72,123.35,132.19,151.87,161.5,177.34,186.92,201.35,216.09,230.12,245.47,255.85,273.45,285.91,303.99,315.98,325.48,343.01,360.05,373.17,381.7,398.41,412.66,423.66,443.67,450.39,468.86,483.93,499.91,511.59,529.34,541.35,550.28,568.31,584.7,592.33,615.74,622.45,639.1,651.41,668.08,679.75,692.94,708.83,720.98,734.42,747.83,762.27,778.74,790.97,806.99,820.03,831.55,844.23),
                   group=factor(rep(c('A','B'),c(81,110))));

plot(y~x,updown,subset=x<500,col=group);

Nuage de points

Vous pourriez d'abord essayer une transformation de Box-Cox et voir si cela avait un sens mécanique, et à défaut, vous pourriez adapter un modèle des moindres carrés non linéaires avec une fonction de lien logistique ou asymptotique.

Alors, pourquoi devriez-vous abandonner complètement les modèles paramétriques et vous rabattre sur une méthode de boîte noire comme les splines lorsque vous découvrez que l'ensemble de données complet ressemble à ceci ...

plot(y~x,updown,col=group);

Mes questions sont:

  • Quels termes dois-je rechercher pour trouver des fonctions de lien qui représentent cette classe de relations fonctionnelles?

ou

  • Que dois-je lire et / ou rechercher pour m'apprendre à concevoir des fonctions de lien avec cette classe de relations fonctionnelles ou à étendre celles existantes qui ne sont actuellement que pour des réponses monotones?

ou

  • Heck, même quelle balise StackExchange est la plus appropriée pour ce type de question!
f1r3br4nd
la source
4
Je n'ai aucune idée de ce que tu demandes. Vous voulez adapter une fonction non monotone de ... quel est exactement votre problème avec la régression polynomiale ou la régression sinusoïdale ?? Aussi ... "fonction de lien" ... vous continuez à utiliser ce mot ... Je ne pense pas que cela signifie ce que vous pensez qu'il signifie. x
Jake Westfall
5
(1) Votre Rcode contient des erreurs de syntaxe: groupne doit pas être cité. (2) L'intrigue est magnifique: les points rouges présentent une relation linéaire tandis que les noirs peuvent être ajustés de plusieurs façons, y compris une régression linéaire par morceaux (obtenue avec un modèle de point de changement) et peut-être même comme exponentielle. Cependant, je ne les recommande pas , car les choix de modélisation doivent être éclairés par une compréhension de ce qui a produit les données et motivés par les théories dans les disciplines pertinentes. Ils pourraient être un meilleur début pour votre recherche.
whuber
1
@whuber merci! Correction du code. Concernant la motivation théorique: d'où viennent-elles en premier lieu? Mes collaborateurs de laboratoire seront heureux de dichotomiser les variables prédictives et d'effectuer des tests t sur celles-ci. Il m'incombe donc de trouver un moyen d'arrêter de gaspiller des données en trouvant une relation mathématique qui capture la transition de "y corrèle positivement avec x" à "y a peu de réponse à x" à "y corrèle négativement avec x". À défaut, je devrai récapituler ce que, par exemple, Michaelis et Menten ont fait lorsqu'ils ont trouvé une relation entre l'enzyme, le substrat et le produit.
f1r3br4nd
1
Est-ce que les points où ces choses se «tordent» sont connus à l'avance?
Glen_b -Reinstate Monica
3
+1 pour le titre provocateur et un suivi qui a du sens
Stumpy Joe Pete

Réponses:

45

Les remarques dans la question sur les fonctions de liaison et la monotonie sont un hareng rouge. Sous-jacents, il semble y avoir une hypothèse implicite selon laquelle un modèle linéaire généralisé (GLM), en exprimant l'attente d'une réponse tant que fonction monotone d'une combinaison linéaire de variables explicatives , n'est pas suffisamment flexible pour tenir compte des non réponses -monotoniques. Ce n'est tout simplement pas le cas.f X β XYfXβX


Peut-être qu'un exemple concret éclairera ce point. Dans une étude de 1948 (publiée à titre posthume en 1977 et jamais évaluée par des pairs), J. Tolkien a rapporté les résultats d'une expérience d'arrosage de plantes dans laquelle 13 groupes de 24 tournesols ( Helianthus Gondorensis ) ont reçu des quantités contrôlées d'eau à partir de la germination pendant trois mois de croissance. Les quantités totales appliquées variaient d'un pouce à 25 pouces par incréments de deux pouces.

Figure 1

Il y a une réponse positive claire à l'arrosage et une forte réponse négative à l'arrosage excessif. Des travaux antérieurs, basés sur des modèles cinétiques hypothétiques de transport ionique, avaient émis l'hypothèse que deux mécanismes concurrents pourraient expliquer ce comportement: l'un a entraîné une réponse linéaire à de petites quantités d'eau (mesurées dans les logarithmes de survie), tandis que l'autre- -un facteur inhibiteur - a agi de façon exponentielle (ce qui est un effet fortement non linéaire). Avec de grandes quantités d'eau, le facteur inhibiteur submergerait les effets positifs de l'eau et augmenterait sensiblement la mortalité.

Soit le taux d'inhibition (inconnu) (par unité de quantité d'eau). Ce modèle affirme que le nombre de survivants dans un groupe de taille recevant pouces d'eau devrait avoir un distribution, où est la fonction de lien convertissant les cotes du journal en une probabilité. Il s'agit d'un GLM binomial. En tant que tel, bien qu'il soit manifestement non linéaire en , étant donné toute valeur de il est linéaire dans ses paramètres , etY n x Binôme ( n , f ( β 0 + β 1 x - β 2 exp ( κ x ) ) ) f x κ β 0 β 1 β 2 f - 1 ( E [ Y ] ) x 1 β 0 x β 1 - exp ( κ x ) β 2κYnx

Binomial(n,f(β0+β1xβ2exp(κx)))
fxκβ0β1β2. La "linéarité" dans le paramètre GLM doit être comprise dans le sens où est une combinaison linéaire de ces paramètres dont les coefficients sont connus pour chaque . Et ils sont: ils sont égaux à (le coefficient de ), lui-même (le coefficient de ) et (le coefficient de ).f1(E[Y])x1β0xβ1exp(κx)β2

Ce modèle - bien qu'il soit quelque peu nouveau et pas complètement linéaire dans ses paramètres - peut être adapté à l'aide d'un logiciel standard en maximisant la probabilité d'arbitraire et en sélectionnant le pour lequel ce maximum est le plus grand. Voici le code pour le faire, en commençant par les données:κκκR

water <- seq(1, 25, length.out=13)
n.survived <- c(0, 3, 4, 12, 18, 21, 23, 24, 22, 23, 18, 3, 2)
pop <- 24
counts <- cbind(n.survived, n.died=pop-n.survived)
f <- function(k) {
  fit <- glm(counts ~ water + I(-exp(water * k)), family=binomial)
  list(AIC=AIC(fit), fit=fit)
}
k.est <- optim(0.1, function(k) f(k)$AIC, method="Brent", lower=0, upper=1)$par
fit <- f(k.est)$fit

Il n'y a pas de difficultés techniques; le calcul ne prend que 1/30 seconde.

Figure 2

La courbe bleue est l'espérance ajustée de la réponse, E[Y] .

E[Y]xR

x.0 <- seq(min(water), max(water), length.out=100)
p.0 <- cbind(rep(1, length(x.0)), x.0, -exp(k.est * x.0))
logistic <- function(x) 1 - 1/(1 + exp(x))
predicted <- pop * logistic(p.0 %*% coef(fit))

plot(water, n.survived / pop, main="Data and Fit",
     xlab="Total water (inches)", 
     ylab="Proportion surviving at 3 months")
lines(x.0, predicted / pop, col="#a0a0ff", lwd=2)

Les réponses aux questions sont:

Quels termes dois-je rechercher pour trouver des fonctions de lien qui représentent cette classe de relations fonctionnelles?

Aucun : ce n'est pas le but de la fonction de liaison.

Que dois-je ... rechercher pour ... étendre les [fonctions de liaison] existantes qui ne sont actuellement réservées qu'aux réponses monotones?

Rien : cela est basé sur une mauvaise compréhension de la façon dont les réponses sont modélisées.

Évidemment, il faut d'abord se concentrer sur les variables explicatives à utiliser ou à construire lors de la construction d'un modèle de régression. Comme suggéré dans cet exemple, recherchez les conseils de l'expérience et de la théorie passées.

whuber
la source
réponse géniale! Ces données sont-elles tolkiennes du roman?
Cam.Davidson.Pilon
1
@Cam Les données ne sont pas entrées dans la coupe finale :-). (Le contexte est plutôt
ironique
1
κ
5
κκχ2(1)
1
@zipzapboing L'exemple que je donne ici est spécial car il a été informé par une théorie sous-jacente. Lorsque de telles informations sont disponibles, elles peuvent être un puissant guide pour sélectionner un modèle. Dans de nombreux cas, cependant, il n'y a pas de telles informations, ou on espère seulement que la réponse attendue pourrait varier de façon monotone avec les régresseurs. La raison la plus fondamentale à laquelle on pourrait faire allusion est peut-être l'espoir que la réponse varie différemment selon les régresseurs et que, pour la plage de régresseurs dans les données, le changement de dérivée est faible: une réponse linéaire se rapprocherait aussi bien de cela.
whuber
9

Regarde coupablement l'usine mourante sur son bureau ... apparemment pas

Dans les commentaires, @whuber dit que "les choix de modélisation doivent être éclairés par une compréhension de ce qui a produit les données et motivés par les théories dans les disciplines pertinentes", à laquelle vous avez demandé comment procéder.

La cinétique de Michaelis et Menten est en fait un exemple assez utile. Ces équations peuvent être dérivées en partant de certaines hypothèses (par exemple, le substrat est en équilibre avec son complexe, l'enzyme n'est pas consommée) et certains principes connus (la loi de l'action de masse). Murray's Mathematical Biology: An Introduction parcourt la dérivation du chapitre 6 (je parierais que beaucoup d'autres livres le font aussi!).

Plus généralement, il permet de constituer un "répertoire" de modèles et d'hypothèses. Je suis sûr que votre domaine a des modèles éprouvés et couramment acceptés. Par exemple, si quelque chose se charge ou se décharge, j'atteindrais une exponentielle pour modéliser sa tension en fonction du temps. Inversement, si je vois une forme exponentielle dans un tracé tension-temps, ma première supposition serait que quelque chose dans le circuit se décharge capacitivement et, si je ne savais pas ce que c'était, j'essaierais de le trouver. Idéalement, la théorie peut à la fois vous aider à construire le modèle et suggérer de nouvelles expériences.

y=k(x+h)2CO2 capture à partir de moins de transpiration?) et les inondations (bactéries mangeant les racines?) pourraient suggérer une forme spécifique pour chaque pièce.

Matt Krause
la source
8

J'ai une réponse plutôt informelle du point de vue de quelqu'un qui a passé la moitié de sa vie scientifique au banc et l'autre moitié à l'ordinateur, à jouer avec les statistiques. J'ai essayé de faire un commentaire, mais c'était trop long.

Vous voyez, si j'étais un scientifique observant le type de résultats que vous obtenez, je serais ravi. Les diverses relations monotones sont ennuyeuses et difficiles à distinguer. Cependant, le type de relation que vous nous montrez suggère un effet très particulier. Cela nous donne un merveilleux terrain de jeu pour le théoricien pour proposer des hypothèses sur ce qu'est la relation, comment elle change aux extrêmes. Cela donne un excellent terrain de jeu au scientifique du banc pour comprendre ce qui se passe et expérimenter largement sur les conditions.

Dans un sens, je préfère avoir le cas que vous montrez et ne pas savoir comment adapter un modèle simple (mais être en mesure d'élaborer une nouvelle hypothèse) que d'avoir une relation simple, facile à modéliser mais plus difficile à étudier mécaniquement. Cependant, je n'ai pas encore rencontré un tel cas dans ma pratique.

Enfin, il y a encore une considération. Si vous cherchez un test qui montre que le noir est différent du rouge (dans vos données) - en tant qu'ancien scientifique de laboratoire, je dis pourquoi même s'embêter? C'est assez clair d'après la figure.

janvier
la source
5

Pour des données comme ça, je considérerais probablement au moins les splines linéaires.

Vous pouvez les faire assez facilement en lm ou glm.

Si vous adoptez une telle approche, votre problème sera de choisir le nombre de nœuds et les emplacements des nœuds; une solution pourrait être de considérer un bon nombre d'emplacements possibles et d'utiliser quelque chose comme le lasso ou d'autres méthodes de régularisation et de sélection pour identifier un petit ensemble; vous devrez cependant prendre en compte l'effet d'une telle sélection dans l'inférence.

Glen_b -Reinstate Monica
la source
Mais la régression spline ne dit-elle pas fondamentalement "il existe une fonction inconnue décrivant la forme de la réponse et nous ne testerons que des hypothèses sur la façon dont les autres variables déplacent cette courbe vers le haut / bas ou l'inclinent"? Et si un traitement modifie la forme elle-même - comment interpréter un tel terme d'interaction s'il est significatif?
f1r3br4nd
2
Quelle est l'alternative générale? Même pour le cas général, il existe une variété d'approches où vous pouvez comparer l'ajustement en supposant des fonctions non paramétriques identiques par rapport à des fonctions distinctes. Les modèles additifs et les modèles additifs généralisés peuvent traiter de telles comparaisons.
Glen_b -Reinstate Monica
À titre d'exemple d'un cas plus général que celui que vous discutez (avec des références discutant d'une variété d'autres approches), si vous pouvez vous en procurer, jetez un œil à cet article J.Roca-Pardiñas et al (2006) "Bootstrap-based méthodes pour tester les interactions facteur par courbe dans des modèles additifs généralisés: évaluer l'activité neuronale du cortex préfrontal liée à la prise de décision ", Statistics in Medicine , 30 juil. 25 (14): 2483-501. Dans cet article, ils utilisent le bootstrap (et le binning pour réduire la charge de calcul), mais il y a d'autres approches mentionnées ici.
Glen_b -Reinstate Monica
Une référence plus basique et plus ancienne serait quelque chose comme Hastie et Tibshirani (1990), Generalized Additive Models (par exemple, voir p265). De plus, jetez un oeil ici , plus précisément, la dernière équation de la diapositive 34. Autour de là , il explique aussi comment adapter un tel modèle à l' aide gamdans le package R mgcv.
Glen_b -Reinstate Monica
2

Je n'ai pas eu le temps de lire l'intégralité de votre article, mais il semble que votre principale préoccupation soit que les formes fonctionnelles des réponses puissent changer avec les traitements. Il existe des techniques pour y faire face, mais elles sont gourmandes en données.
Pour votre exemple spécifique:

G est la croissance W est l'eau T est le traitement

library(mgcv)
mod = gam(G~T+s(W,by=T))
plot(mod,pages=1,all=TRUE)
?gam

La dernière décennie a vu une tonne de recherches sur la régression semi-paramétrique, et ces boeufs sur les formes fonctionnelles deviennent de plus en plus gérables. Mais à la fin de la journée, les statistiques jouent avec les chiffres et ne sont utiles que dans la mesure où elles renforcent l'intuition des phénomènes observés. Cela nécessite à son tour de comprendre la façon dont les chiffres sont joués. Le ton de votre message indique une volonté de jeter le bébé avec l'eau du bain.

utilisateur_générique
la source