La régression logistique est comme la régression classique mais permet cette fois de regarder le lien entre une variable qualitative et une variable quantitative. Prenons par exemple nos données sur les rapports entre l'état d'une tumeur sur la vessie et l'expression des gènes que nous avons choisis. Il s'agit cette fois d'expliquer comment l'expression d'un gène influe sur la probabilité que la tumeur soit dans un des états que le chirurgien a observé dans son échantillon.
Supposons dans un premier temps que l'on se limite au cas "invasif contre superficiel". On veut donc savoir si l'expression d'un gène intervient dans la probabilité que la tumeur soit superficielle ou invasive. Comme il n'y a que deux niveaux sur le facteur "Cl.Clinique", il suffit de se fixer un probabilité p pour que l'état soit invasif. Ainsi, l'état sera superficiel avec probabilité 1-p.
Comment relier l'expression X du gène considéré à p ? il faudrait que la relation fasse intervenir quelque chose du type aX+b comme pour la regression linéaire, afin que nous puissons interpréter le coefficient a et tester s'il est nul ou pas pour affirmer ou non que le gène étudié est influent ou non.
Les gens ont réalisé un consensus sur une relation du type
p=exponentielle (aX+b)/(1+exponentielle(aX+b))
Dans le cas de la régression linéaire classique, les coefficients a et b étaient estimés par "moindres carrés". Une methode pour estimer les coefficients a et b dans notre cas logistique classique est une généralisation appelée méthode du maximum de vraisemblance. Elle correspond a chercher a et b qui rendent les plus probables possible, avec le modèle correspondant, les données observées. Une procédure toute faite est donnée sous R: la procedure "glm"
il suffit de faire
> choupette <- glm(vessie$Cl.clinique~vessie$GEC1, family=binomial(link="logit"))
>
> summary(choupette)
Call:
glm(formula = vessie$Cl.clinique ~ vessie$GEC1, family = binomial(link = "logit"),
na.action = na.pass)
Deviance Residuals:
Min 1Q Median 3Q Max
-1.8849 0.6220 0.6748 0.7325 1.6816
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 1.62859 0.38280 4.254 2.10e-05 ***
vessie$GEC1 -0.16334 0.09974 -1.638 0.101
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 110.22 on 99 degrees of freedom
Residual deviance: 107.34 on 98 degrees of freedom
AIC: 111.34
Number of Fisher Scoring iterations: 4
Comment traduire toutes informations obscures en quelque chose de pertinent ? Comme
pour la régression linéaire classique, il faut trouver dans ce galimatia les valeurs
estimées pour a et b et la p-value pour le test qui permet de savoir si H0: a=0 ou
H1: a différent de 0 est vraie.
On voit que les coefficients sont donnés comme pour la régression linéaire par
> choupette$coefficients
(Intercept) vessie$GEC1
1.6285933 -0.1633434
Le test de nullité pour a, si important afin de savoir si le gène en question est influent par rapport au stade de la tumeur est donné dans le summary sous le signen Pr(>|z|) et est 0.101 pour le coefficient
a comme vous pouvez le vérifier. On conclue que a n'est pas égal à 0 et que le gène GGEC1 intervient !
Il est important dans certaine situations comme dans la situation présente de prendre en compte le fait que tous les gènes influent en même temps sur la réponse, ici le stade d'evolution de la tumeur. Soient par exemple Y et Z les variables d'expressio de deux autres gènes d'intérêt. On peut faire un modèle logistique un peu plus compliqué qui prend en compte l'influence simultannnée de ces trois gènes. Le modèle est comme suit:
p=exponentielle (a_1 X+a_2 Y+a_3 Z+b)/(1+exponentielle(a_1 X+a_2 Y+a_3 Z+b))
> choupette <- glm(vessie$Cl.clinique~vessie$GEC1+vessie$UBE2C+vessie$GEC3, family=binomial(link="logit"))
Les coefficients sont donnés par la requêtte suivante:
> choupette$coefficients
(Intercept) vessie$GEC1 vessie$UBE2C vessie$GEC3
2.41683217 -0.07707264 -6.39926485 0.17864624
Toutes les informations sur le résultat de cette estimation sont obtenu en tapant "summary" comme d'habitude:
> summary(choupette)
Call:
glm(formula = vessie$Cl.clinique ~ vessie$GEC1 + vessie$UBE2C +
vessie$GEC3, family = binomial(link = "logit"), na.action = na.pass)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.1363 0.2101 0.4284 0.5461 1.9900
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 2.41683 0.54794 4.411 1.03e-05 ***
vessie$GEC1 -0.07707 0.10643 -0.724 0.469
vessie$UBE2C -6.39926 1.54328 -4.147 3.38e-05 ***
vessie$GEC3 0.17865 0.12980 1.376 0.169
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 110.216 on 99 degrees of freedom
Residual deviance: 82.387 on 96 degrees of freedom
AIC: 90.387
Number of Fisher Scoring iterations: 5
On interprète ces résultats facilement maintenant qu'on a l'habitude: l'estimation de a_1 donne -0.07707, celle de a_2 donne -6.39926 et pour a_3 on obtient 0.17865. Les p-values associées sont
respectivement 0.469, 3.38e-05 et 0.169, ce qui donne que UBE2C est très significativement influent alors que les GEC le sont beaucoup moins !
Une conclusion est qu'il ne faut jamais hésiter à mettre beaucoup de variables explicatives dans le même panier pour expliquer une réponse. Ainsi, chaque variable explique sa part d'effet sur la sortie observée. De part les corrélations entre les variables explicatives, on va de plus obtenir des choses différentes si on explique la sortie observée selon chaque co-variable l'une après l'autre ou toutes ensemble ...
Le hic est que si on n'a que 100 patients et 100 gènes par exemple, ce la fait 100 coefficients de régression a_1, a_2, a_3 etc ... à estimer et 1 intersept b, donc 101 inconnues pour 100 observations ! on est souvent dans des situations encore plus abracadabrantes ou le nombre de paramètres à estimer est 10 fois plus grand que le nombre de patients observés ... Pour s'en sortir, une seule issue: profiter de la sparsité (c'est a dire: en vrai peu de gènes sont influent pour une maladie données et donc la plupart des a_i sont nuls sans qu'on sache à l'avance lesquels), et utiliser des techniques comme le LASSO ou le LASSO logistique. Pour l'instant, ce sont surtout des voies de recherche mais cela devrait pouvoir être utilisable en pratique trềs bientôt !