La régression des moindres carrés partiels

add2.png

\bullet Présentation:

La régression des moindres carrés partiels a été élaborée en 1983 par Herman Wold et Svante Wold et est souvent retrouvée sous l’appellation de régression PLS (Partial Least Squarre). Elle permet de discriminer une variable réponse Y continue, binaire ou polychotomique (K \geq 2 classes) à partir d’une matrice de P variables explicatives \mathbf{X} = (X ^1, \cdots, X ^P) continues, qualitatives ou mixtes.

L’outil a été développé afin de contourner les deux principales problématiques rencontrées par la régression linéaire et logistique:

– La multicolinéarité qui induit en erreur l’algorithme de construction de la règle décisionnelle en dispersant l’information d’une variable sur celles avec lesquelles elle est corrélée.

– La présence de beaucoup plus de variables que d’individus (le fameux problème de type P >>>>> N, amenant un biais dans l’estimation due à une trop grande variance.

L’objectif de la régression PLS étant de combiner les caractéristiques de l’analyse en composantes principales (ACP) et celles de la régression afin de réaliser un compromis entre maximisation de la variance expliquée par les variables \mathbf{X} et maximisation de leur corrélation avec Y.

\bullet La régression PLS:

Hypothèse préliminaire: Variables explicatives \mathbf{X} continues et qualitatives. Variable réponse Y continue, binaire ou polychotomique K \geq 2.

L’Algorithme:

De manière générale, les données continues doivent être centrées-réduites. L’Algorithme présenté ci-dessous, inspiré de l’algorithme NIPALS, implémente le calcul étape par étape des composantes D_h. En fonction de la version de la régression (logistique ou linéaire) PLS, la méthode de modélisation varie. A noter que la régression PLS nécessite le paramètrage du nombre H de composantes à calculer.

– Étape 1: Calcul de D_1

♦ Calcul du coefficient w_1 ^j, associé à la variable X ^j, par régression logistique ou linéaire univariée de Y, \forall j \in \lbrace 1, \cdots, P \rbrace.

♦ Normalisation du vecteur des coefficients w_1 en w_{1,N} = \frac{w_1}{\sqrt{\sum_{j = 1} ^P (w_1 ^j) ^2}}.

♦ Calcul de la composante D_1 = \frac{X \cdot w_{1,N}}{w_{1,N} \cdot w_{1,N} ^t}.

♦ Détermination du vecteur des coefficients w_{1,*} par régression linéaire de D_1 sur \mathbf{X} sous contrainte que \vert w_{1,*} \vert = 1 et sans estimation du paramètre constant.

♦ Calcul de la matrice résiduelle \mathbf{X}_1 via la formule: \mathbf{X}_1 = \mathbf{X} - D_1 \cdot (\frac{D_1 ^t \cdot \mathbf{X}}{D_1 ^t \cdot D_1}).

♦ Pour la version linéaire de la régression PLS exclusivement, nous devons également calculer le vecteur résiduel Y_1 via la formule: Y_1 = Y - D_1 \cdot (\frac{D_1 ^t \cdot Y}{D_1 ^t \cdot D_1}).

– Étape h \in [2, H]: Calcul de D_h

♦ Calcul du coefficient w_h ^j associé à la variable X_{h-1} ^j dans la régression multiple logistique  de Y ou linéaire de Y_{h-1} sur les variables D_1, \ldots, D_{h-1}, X_{h-1} ^j, \forall j \in \lbrace 1, \ldots, P \rbrace.

♦ Normalisation du vecteur des coefficients w_h en w_{h,N} = \frac{w_h}{\sqrt{\sum_{j = 1} ^P (w_h ^j) ^2}}.

♦ Calcul la composante D_h = \frac{\mathbf{X}_{h-1} \cdot w_{h,N}}{w_{h,N} \cdot w_{h,N} ^t}.

♦ Détermination du vecteur des coefficients w_{h,*} par régression linéaire de D_h sur \mathbf{X} sous contrainte que \vert w_{1,*} \vert = 1 et sans estimation du paramètre constant.

♦ Calcul de la matrice résiduelle \mathbf{X}_h via la formule: \mathbf{X}_h = \mathbf{X}_{h-1} - D_h \cdot (\frac{D_h ^t \cdot \mathbf{X}_{h-1}}{D_h ^t \cdot D_h}).

♦ Pour la version linéaire de la régression PLS exclusivement, nous devons également calculer le vecteur résiduel Y_h via la formule: Y_h = Y_{h-1} - D_h \cdot (\frac{D_h ^t \cdot Y_{h-1}}{D_h ^t \cdot D_h}).

– Étape finale: Calcul des coefficients de la règle décisionnelle

Nous nous retrouvons à la fin de la H étape de l’algorithme avec nos composantes D_1, \cdots, D_H et nos vecteurs de coefficients w_{1,*}, \cdots, w_{H,*}.

♦ Dans un premier temps nous estimons les coefficients \alpha_0, \cdots, \alpha_h par régression logistique ou linéaire standard de Y sur D_1, \cdots, D_h.

Nous avons alors le coefficient \beta_0 = \alpha_0 et sera conservé comme le coefficient constant du modèle final.

♦ Enfin, la détermination des coefficients finaux \beta_1, \dots, \beta_P se fait en développant l’équation suivante:

(\beta_1, \dots, \beta_P) = \alpha_1 \dot w_{1,*} + \cdots + \alpha_h \cdot w_{h,*}

La règle décisionnelle est la même que celle de la régression logistique ou linéaire, en fonction de la version choisie, et s’exprime en fonction du vecteur des coefficients (\beta_0, \cdots, \beta_P).

Caractéristiques:

– Le calcul des composantes D_h, h \in [1,H] est analogue à celui des composantes principales de l’ACP à ceci prés que dans le cas de la régression PLS elles sont déterminées sous contrainte de maximiser la covariance entre les variables explicatives et la variable réponse.

– La problématique de la multicolinéarité est traitée en deux temps: lors de la phase de régression univariée qui revient à calculer les coefficients « individuellement » de chaque variable et donc à retenir l’information la plus brut qui soit. Et lors de la régression sur les composantes D_h, h \in [1,H] qui synthétisent l’information maximisée par projection orthogonale et la fusionne.

– La gestion du problème de type P >>>>>> N est traitée lorsque de la phase de régression univariée et limite ainsi les risques de sur-apprentissage du modèle.

– En pratique, le nombre de composantes Dh nécessaires et suffisantes est rarement supérieur à trois ou quatre. Au delà les performances n’évoluent quasiment plus. Ce dernier est choisi soit par validation-croisée soit par maximisation des performances.

– L’Algorithme PLS est très rapide étant donné qu’il consiste en une succession de régressions d’où son efficacité sur de grandes bases de données et ne nécessite aucune opération matricielle complexe.

– La régression PLS tend à devenir la modélisation dominante et reste généralement plus efficace que la régression logistique, la régression linéaire, la régression sur composantes principales et la régression Ridge.

Dans le cas de la régression logistique, les coefficients \beta donne la valeur des odds ratios lorsque nous leur appliquons la fonction exponentielle. Le calcul des intervalles de confiance pour le rejet de l’hypothèse \beta_j = 0 se fait par bootstrap.

\bullet Annexe théorique:

Nous présentons ici la démonstration de la recherche de la première composante D_1.

Nous avons vu que la composante D_1 peut s’écrire de la manière suivante,

D_1 = \mathbf{X} w_{1,*} = \mathbf{X} w_1

Elle est alors obtenue en maximisant le critère de Tucker,

cov ^2 (Y, \mathbf{X} w_1) = r ^2 (Y, \mathbf{X} w_1)) Var( \mathbf{X} w_1) Var(Y)

, sous contrainte que ||w_1|| = 1.

L’idée est d’essayer de maximiser parallèlement la variance expliquée par D_1 et la corrélation entre D_1, Y tel que,

<D_1, Y> = <\mathbf{X} w_1, Y> = || \mathbf{X} w_1 || \cdot ||Y|| \cdot cor(\mathbf{X} w_1, Y)

Le vecteur des coefficients w_1 est déterminé à l’aide des multiplicateurs Lagrangien \lambda \in R_+,

L(w_1, \lambda) = cov(Y,\mathbf{X} w_1) - \lambda(w_1 ^t w_1 - 1) = w_1 ^t \mathbf{X} ^t Y - \lambda(w_1 ^t w_1 - 1)

Il faut ensuite annuler les dérivées,

\frac{\partial L}{\partial \lambda} = - (w_1 ^t w_1 - 1) = 0 \Rightarrow w_1 ^t w_1 = 1

\frac{\partial L}{\partial w_1} = \mathbf{X} ^t Y - 2 \lambda w_1 = 0 \Rightarrow \mathbf{X} ^t Y = 2 \lambda w_1

Nous avons alors,

w_1 ^t \mathbf{X} ^t Y = w_1 ^t 2 \lambda w_1 = 2 \lambda w_1 ^t w_1 = 2 \lambda

Maintenant, soit \theta \in R tq,

\theta = 2 \lambda = Y ^t \mathbf{X} w_1 \Rightarrow (\mathbf{X} ^t Y)(Y ^t \mathbf{X} w_1 = (\theta w_1) \theta = \theta ^2 w_1

Nous en déduisons que w_1 est le vecteur propre de \mathbf{X} ^t Y Y ^t \mathbf{X} associé à la valeur propre \theta ^2

Nous pouvons déduire de ce résultat que,

\lambda_1 = <\mathbf{X} w_1, Y> ^2 = Y ^t \mathbf{X} \mathbf{X} ^t Y

Par conséquent,

\mathbf{X} ^t Y Y ^t \mathbf{X} w_1 = (\mathbf{X} ^t Y) Y ^t \mathbf{X} w_1 = \lambda_1 w_1 = Y ^t \mathbf{X} \mathbf{X} ^t Y w_1 = w_1 Y ^t \mathbf{X} \mathbf{X} ^t Y

\Rightarrow w_1 = \mathbf{X} ^t Y

La division par la norme est alors obligatoire sous contrainte fixée sur w_1.

\bullet Exemple:

Soit l’échantillon suivant,

add

La variable « Statut » représente la réponse et les variables X ^1, X ^2 nos variables explicatives.

Le travail préliminaire à effectuer étant de transformer la matrice \mathbf{X} = (X ^1, X ^2) en un tableau disjonctif complet pour lequel nous supprimons les colonnes correspondantes à la modalité « 1m » pour nos deux variables explicatives. Cette modalité faisant office de modalité de référence.

Nous obtenons alors une matrice \mathbf{X} de taille 4 \times 20 étant donné que nous avons désormais quatre variables: X_{1m} ^1, X_{2m} ^1, X_{1m} ^2, X_{2m} ^2.

– Etape 1: calcul de la composante D_1.

Pour débuter, nous procédons donc à quatre régression logistique univariée de Y sur X_{1m} ^1, puis sur X_{2m} ^1, puis sur X_{1m} ^2 et enfin sur X_{2m} ^2. Nous obtenons  ainsi,

w_1 = (-0.4418, 0.4418, -0.4418 , 2.6027)

Nous calculons maintenant le coefficient de normalisation,

\sqrt{(-0.4418)^2 + 0.4418^2 + (-0.4418)^2 + 2.6027^2} = 2.71286

Nous obtenons alors le vecteur des coefficients normalisés,

w_{1,N} = \frac{w_1}{2.71286} = (-0.1628540, 0.1628540, -0.1628540, 0.9593934)

Désormais, nous pouvons calculer la composante:

D_1 = \frac{X \cdot w_{1,N}}{w_{1,N} \cdot w_{1,N} ^t}

 = \frac{\begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 1 & 0 & 1 & 0 \\ 0 & 1 & 1 & 0 \\ 1 & 0 & 1 & 0 \\ 1 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \\ 0 & 1 & 0 & 0 \\ \vdots & \vdots & \vdots & \vdots \end{pmatrix} \times (-0.1628540, 0.1628540, -0.1628540, 0.9593934)}{(-0.1628540, 0.1628540, -0.1628540, 0.9593934) \times \begin{pmatrix} -0.1628540 \\ 0.1628540 \\ -0.1628540 \\ 0.9593934 \end{pmatrix}}

 = (0, 0, -0.1628540, 0.1628540, 0.1628540, -0.3257079, 0, -0.3257079, -0.3257079, 0.9593934, 0.1628540,

-0.3257079, 0, -0.1628540, 0.7965395, 0.9593934, 0.9593934, 1.1222474, 0.7965395, 1.1222474) \times \frac{1}{1}

 = (0, 0, -0.1628540, 0.1628540, 0.1628540, -0.3257079, 0, -0.3257079, -0.3257079, 0.9593934,

0.1628540, -0.3257079, 0, -0.1628540, 0.7965395, 0.9593934, 0.9593934, 1.1222474, 0.7965395, 1.1222474)

Ensuite, nous devons calculer le vecteur des coefficients w_{1,*}. Pour cela nous appliquer une régression linéaire de D_1 sur X et obtenons le vecteur des coefficients suivant,

w_{1,*} = (-0.1628657,0.1628657,-0.1628657,0.9593875)

Enfin, nous concluons l’étape 1 en déterminant la matrice résiduelle \mathbf{X}_1. Nous appliquons donc la formule:

\mathbf{X}_1 = \mathbf{X} - D_1 \cdot \frac{D_1 ^t \cdot \mathbf{X}}{D_1 ^t \cdot D_1}

Nous calculons dans un premier temps,

\frac{D_1 ^t \cdot X \cdot}{D_1 ^t \cdot D_1} = \frac{(0.1273932, 2.733057, -1.465686, 6.715754)}{7.106086} = (0.01792734, 0.38460788, -0.20625780, 0.94507075)

Et,

D_1 \cdot (0.01792734, 0.38460788, -0.20625780, 0.94507075) = \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 1.002919539 & 0.06263492 & -0.0335899 & 0.15390852 \\ -0.002919539 & 0.93736508 & 0.0335899 & -0.15390852 \\ -0.002919539 & 0.93736508 & 0.0335899 & -0.15390852 \\ 1.005839077 & 0.12526984 & 0.9328202 & 0.30781704 \\ 0 & 1 & 1 & 0 \\ 1.005839077 & 0.12526984 & 0.9328202 & 0.30781704 \\ 1.005839077 & 0.12526984 & 0.9328202 & 0.30781704 \\ -0.017199374 & -0.36899028 & 0.1978824 & 0.09330533 \\ -0.002919539 & 0.93736508 & 0.0355899 & -0.15390852 \\ \vdots & \vdots & \vdots & \vdots \end{pmatrix}

Il nous reste plus qu’à soustraire ce résultat à \mathbf{X} selon la formule de calcul de la matrice résiduelle \mathbf{X}_1. Nous trouvons alors,

\mathbf{X}_1 = \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 1.002919539 & 0.06263492 & -0.0335899 & 0.15390852 \\ -0.002919539 & 0.93736508 & 0.0335899 & -0.15390852 \\ 0.000000000 & 1.00000000 & 1.0000000 & 0.00000000 \\ 0.002919539 & 0.06263492 & 0.9664101 & 0.15390852 \\ -0.020118913 & 0.56837480 & 0.2314723 & -0.06060319 \\ -0.017199374 & -0.36899028 & 0.1978824 & 0.09330533 \\ -0.002919539 & 0.93736508 & 0.0335899 & -0.15390852 \\ 0.000000000 & 0.00000000 & 0.0000000 & 0.00000000 \\ 0.000000000 & 0.00000000 & 0.0000000 & 0.00000000 \\ 1.005839077 & 0.12526984 & 0.9328202 & 0.30781704 \\ \vdots & \vdots & \vdots & \vdots \end{pmatrix}

– Etape 2: calcul de la composante D_2.

L’étape 2 est similaire à l’étape 1 à l’exception de la détermination du vecteur des coefficients w_2 où nous ne faisons plus des régressions logistiques univariées mais des régressions logistiques multivariées de (D_1, X_1 ^j), \forall j \in [1, 4].

Nous obtenons alors pour la régression logistique de Y sur (D_1, X_{1m} ^1), le coefficient associé à X_{1m} ^1: w_2 ^1 = 0.4234377.

Pour la régression logistique de Y sur (D_1, X_{2m} ^1), le coefficient retenu est: w_2 ^2 = 0.2025917.

Pour la régression logistique de Y sur (D_1, X_{1m} ^2), le coefficient retenu est: w_2 ^3 = 2.4604802.

Pour la régression logistique de Y sur (D_1, X_{2m} ^2), le coefficient retenu est: w_2 ^4 = 2.3654113.

Ce qui nous donne,

w_2 = (0.4234377, 0.2025917, 2.4604802, 2.3654113)

En normalisant le vecteur w_2 nous avons,

w_{2,N} = (0.12290620, 0.05880388, 0.71417409, 0.68657959)

Et la composante D_2,

D_2 = (0, 0, 0.20863288, -0.02692281, -0.02692281, 1.00853366, 0.77297797, 1.00853366, 1.00853366, 0.18159222,

-0.02692281, 1.0085366, 0.77297797, 0.79990078, 0.39022511, 0.18159222, 0.18159222, 0.15466942, 0.39022511, 0.15466942)

Il nous reste plus qu’à déterminer les coefficients w_{2,*}. La régression linéaire de D_2 sur \mathbf{X} nous permet d’obtenir,

w_{2,*} = (0.20863288, -0.02692281, 0.79990078, 0.18159222)

Étant donné que nous avons choisi de nous restreindre à deux composantes, il est inutile de calculer la matrice résiduelle \mathbf{X}_2.

– Etape finale: calcul des coefficients finaux

Nous procédons à l’estimation des coefficients \alpha via régression logistique de Y sur les composantes D_1 et D_2. Nous obtenons alors,

\alpha_0 = -1.862, \alpha_1 = 3.405 et \alpha_2 = 2.420

Ainsi, nous avons déjà le coefficient \beta_0 = \alpha_0 du modèle final. Il nous reste désormais à déterminer les coefficients \beta_1, \beta_2, \beta_3, \beta_4 finaux,

\begin{pmatrix} \beta_1 \\ \beta_2 \\ \beta_3 \\ \beta_4 \end{pmatrix} = w_{1,*} \times \alpha_1 + w_{2,*} \cdot \alpha_2

= \begin{pmatrix} -0.1628657 \\ 0.1628657 \\ -0.1628657 \\ 0.9593875 \end{pmatrix} \times 3.405 + \begin{pmatrix} 0.20863288 \\ -0.02692281 \\ 0.79990078 \\ 0.18159222 \end{pmatrix} \times 2.420

= \begin{pmatrix} -0.04982836 \\ 0.48947624 \\ 1.38074546 \\ 3.70642121 \\ \end{pmatrix}

Nous obtenons donc les coefficients de la régression logistique PLS suivants:

\begin{tabular}{|l|c|c|c|c|c|} \hline Parametres & Coefficients \\ \hline Intercept & -1.86204910 \\ \hline X1 = 2m & -0.04982836 \\ \hline X1 = 3m & 0.48947624 \\ \hline X2 = 2m & 1.38074546 \\ \hline X2 = 3m & 3.70642121 \\ \hline \end{tabular}

En appliquant la formule de l’équation prédictive logit comme en régression logistique standard, mais pour les coefficients calculés via l’algorithme PLS, nous obtenons la matrice de confusion suivante:

\begin{tabular}{|l|c|c|} \hline & Predit Y1 & Predit Y2 \\ \hline Y1 & 7 & 3 \\ \hline Y2 & 4 & 6 \\ \hline \end{tabular}

\bullet Application informatique:

Procédure SAS:

– Régression linéaire PLS: https://support.sas.com/documentation/cdl/en/statug/63033/HTML/default/viewer.htm#statug_pls_sect006.htm

– Régression logistique PLS:

%macro PLS(DATA = , Y = , n_comp = );

/* La macro PLS suivante permet de procéder à la régression logistique des moindres carrés partiels.
Elle nécessite en entrée: – DATA la base de données composées aussi bien de variables continues que qualitatives (elle fonction pour le mixte des deux également)
– Y la variable réponse
– n_comp le nombre de composantes
Elle retourne la table de données mise à jour avec les prédictions et la classification, la table des coefficients ainsi que le tableau de fréquences entre Y et
la prédiction. */

/* Options nécessaires pour éviter de saturer le log et mettre en évidence les erreurs */
options nonotes spool;

/* Récupération du nombre d’observations */
proc sql noprint;
select count(*) into: n from &DATA.;
quit;

/* Vue sur la présence ou non de variable qualitative */
proc contents data = &DATA. (drop =  &Y.) out = biblio noprint;
run;

data biblioQ;
set biblio;
if type = 2;
run;

/* Sortie de la liste des variables qualitatives et de leur nombre */
proc sql noprint;
select distinct(name) into: list_varQ separated by  »  » from biblioQ;
select count(*) into: nb_varQ from biblioQ;
quit;

/* Si des variables qualitatives sont présentes alors on créé le tableau disjonctif qui leur est associé et on supprime la colonne correspondant à la modalité (la première de la liste par défaut */
%if &nb_varQ. > 0 %then %do;

data Xt;
set &DATA.;
drop &list_varQ.;
run;

/* Pour chaque variable sortie du tableau disjonctif avec leur nom correctement associé */
%do p = 1 %to &nb_varQ.;

data encours;
set &DATA.;
keep %scan(&list_varQ.,&p.);
run;

/* Les warnings générés par la macro proviennent d’ici, ils ne sont pas à prendre en compte */
proc corresp data = encours observed binary;
table %scan(&list_varQ.,&p.);
ods exclude binary Inertias RowCoors RowBest RowContr ColBest BinaryPct ColContr ColCoors ColQualMassIn ColSqCos DF RowQualMassIn RowSqCos SupColCoors SupColQuality SupRowCoors SupRowQuality SupRowSqCos;
ods output binary = disj_encours;
run;

data disj_encours;
set disj_encours;
drop Label;
run;

/* Sortie des noms correctes */
proc contents data = disj_encours out = biblioEC noprint;
run;

proc sql noprint;
select distinct(name) into: list_EC separated by  »  » from biblioEC;
select count(*) into: nb_mod from biblioEC;
quit;

/* Inituté correct des noms de variables */
data disj_encours;
set disj_encours;
%do m = 2 %to &nb_mod.;
rename « %scan(&list_EC.,&m.) »n = « %scan(&list_varQ.,&p.)_%scan(&list_EC.,&m.) »n;
%end;
drop « %scan(&list_EC.,1) »n;
run;

/* Ajout à la table Xt = X0 */
data Xt;
merge Xt disj_encours;
run;

%end;

%end;

/* Si aucune variable qualitative alors Xt = X0 = la table d’origine */
%if &nb_varQ. = 0 %then %do;

data Xt;
set &DATA.;
run;

%end;

/* Sauvegarde de l’état modifié de la base pour reconstitution en fin de programme */
data save;
set Xt;
run;

/* Sortie de la nouvelle liste des variables */
proc contents data = Xt (drop =  &Y.) out = biblio noprint;
run;

proc sql noprint;
select distinct(name) into: list_var separated by  »  » from biblio;
select count(*) into: nb_var from biblio;
quit;

/* Création de la table de la variable réponse pour reconstitution en fin de programme */
data Y;
set &DATA.;
keep &Y.;
run;

/* Sortie des classes de la variable réponse pour définir par défaut (la première) la modalité de référence */
proc freq data = Y noprint;
tables &Y. / out = class;
run;

data _null_;
set class;
if _n_ = 1;
call symput(« ref »,&Y.);
run;

/* Pour chaque composante paramétrée on lance l’algorithme PLS */
%do c = 1 %to &n_comp.;

%put Calcul composante N° &c.;

data W;
length Variable $50.;
run;

/* Etape 1 – Sortie des coefficients univariées par régression logistique */
%do p = 1 %to &nb_var.;

data log_univ;
set Xt;
keep &Y. %scan(&list_var.,&p.);
run;

%let model = %scan(&list_var.,&p.);

/* Si on est à la composante 2, 3, 4, etc, alors il ne s’agit plus de régession univariée mais multivariée avec la variable en cours et les composantes T */
%if &c. > 1 %then %do;

%let nb_T = %eval(&c. – 1);

%do t = 1 %to &nb_T.;

data log_univ;
merge log_univ T&t.;
run;

%let model = &model. T&t.;

%end;

%end;

/* Sortie des coefficients de l’étape 1 */
proc logistic data = log_univ;
model &Y. (ref = « &ref. ») = &model.;
ods exclude EffectPlot ROCCurve LeverageCPlot LeverageDifChisqPlot InfluencePlots LeveragePlots LeverageDifDevPlot LeveragePhatPlot CBarPlot CPlot DevianceResidualPlot DifChisqPlot DifDeviancePlot LeveragePlot PearsonResidualPlot ModelInfo ResponseProfile Type3 Classification Association ConvergenceStatus FitStatistics GlobalTests OddsRatios ParameterEstimates ClassLevelInfo NObs;
ods output ParameterEstimates = coeffs;
run;

/* Pour la première composante on garde le coefficient estimé seul, pour les autres on garde seulement celui associé à la variable en cours et ajustée en fonction des composantes T */
data coeffs;
set coeffs;
if Variable = « %scan(&list_var.,&p.) »;
keep Variable Estimate;
run;

/* Ajout à la table des coefficients W */
data W;
set W coeffs;
run;

%end;

/* Etape 2 – Normalisation des coefficients */
data W;
set W;
if Estimate ne .;
Estimate_sqr = Estimate **2;
run;

proc sql noprint;
select sqrt(sum(Estimate_sqr)) into: somE from W;
quit;

/* Application de la formule */
data Wnorm&c.;
set W;
Estimate = Estimate/&somE.;
drop Estimate_sqr;
run;

/* Etape 3 – Calcul des composantes */
data T&c.;
run;

/* Calcul manuel de Xt.W */
%do i = 1 %to &n.;

data Xti;
set Xt;
if _n_ = &i.;
run;

data calc;
run;

%do p = 1 %to &nb_var.;

data _null_;
set Xti;
call symput(« val1 »,%scan(&list_var.,&p.));
run;

data _null_;
set Wnorm&c.;
if _n_ = &p.;
call symput(« val2 »,Estimate);
run;

data calc;
set calc end = eof;
output;
if eof then do;
terme = &val1. * &val2.;
output;
end;
run;

%end;

proc sql noprint;
select sum(terme) into: composante from calc;
quit;

data T&c.;
set T&c. end = eof;
output;
if eof then do;
composante = &composante.;
output;
end;
run;

%end;

data calcD;
run;

/* Calcul manuel de W.W */
%do p = 1 %to &nb_var.;

data _null_;
set Wnorm&c.;
if _n_ = &p.;
call symput(« val »,Estimate);
run;

data calcD;
set calcD end = eof;
output;
if eof then do;
terme = &val. * &val.;
output;
end;
run;

%end;

proc sql noprint;
select sum(terme) into: Denom from calcD;
quit;

/* Calcul manuel de Xt.T/T.T sachant que l’on a calculé les deux termes indépendemment auparavant */
data T&c.;
set T&c.;
if composante ne .;
composante = composante/&Denom.;
run;

/* Etape 4 – Passage de X_t à X_t+1 */
%if &c. < &n_comp. %then %do;

data prod;
run;

/* Calcul manuel de Xt.T */
%do p = 1 %to &nb_var.;

data Xtp;
set Xt;
keep %scan(&list_var.,&p.);
run;

data calcN;
run;

%do i = 1 %to &n.;

data _null_;
set Xtp;
if _n_ = &i.;
call symput(« val1 »,%scan(&list_var.,&p.));
run;

data _null_;
set T&c.;
if _n_ = &i.;
call symput(« val2 »,composante);
run;

data calcN;
set calcN end = eof;
output;
if eof then do;
terme = &val1. * &val2.;
output;
end;
run;

%end;

proc sql noprint;
select sum(terme) into: Num from calcN;
quit;

data prod;
set prod end = eof;
output;
if eof then do;
calc = &Num.;
output;
end;
run;

%end;

data calcD;
run;

   /* Calcul manuel de T.T */
%do i = 1 %to &n.;

data _null_;
set T&c.;
if _n_ = &i.;
call symput(« val »,composante);
run;

data calcD;
set calcD end = eof;
output;
if eof then do;
terme = &val. * &val.;
output;
end;
run;

%end;

proc sql noprint;
select sum(terme) into: Denom from calcD;
quit;

 /* Calcul manuel de Xt.T/T.T sachant que l’on a calculé les deux termes indépendemment auparavant */
data prod;
set prod;
calc = calc/&Denom.;
if calc ne .;
run;

data prodprod;
run;

 /* Calcul manuel de T.(Xt.T/T.T) sachant que l’on a calculé le terme de droite auparavant */
%do i = 1 %to &n.;

data prodprod_ajout;
run;

data _null_;
set T&c.;
if _n_ = &i.;
call symput(« valT »,composante);
run;

%do p = 1 %to &nb_var.;

data _null_;
set prod;
if _n_ = &p.;
call symput(« valP »,calc);
run;

data prodprod_ajout2;
X%scan(&list_var.,&p.) = &valT. * &valP.;
run;

data prodprod_ajout;
merge prodprod_ajout prodprod_ajout2;
run;

%end;

data prodprod;
set prodprod prodprod_ajout;
run;

%end;

data prodprod;
set prodprod;
if X%scan(&list_var.,1) ne .;
run;

data Xt;
merge Xt prodprod;
run;

 /* Calcul manuel de Xt – T.(Xt.T/T.T) */
%do p = 1 %to &nb_var.;

data Xt;
set Xt;
%scan(&list_var.,&p.)b = %scan(&list_var.,&p.) – X%scan(&list_var.,&p.);
drop %scan(&list_var.,&p.) X%scan(&list_var.,&p.);
run;

data Xt;
set Xt;
rename %scan(&list_var.,&p.)b = %scan(&list_var.,&p.);
run;

%end;

data Xt;
merge Y Xt;
run;

%end;

data T&c.;
set T&c.;
rename composante = T&c.;
run;

 /* Etape 5 – Estimation des coefficients w* par régression linéaire */
data reglin;
merge T&c. save (keep = &list_var.);
run;

proc reg data = reglin;
model T&c. = &list_var. / noint;
ods exclude NObs ANOVA FitStatistics ParameterEstimates;
ods output ParameterEstimates = coeff_wetoile&c.;
run;

data coeff_wetoile&c.;
set coeff_wetoile&c.;
keep Variable Estimate;
run;

%end;

/* Etape 6 – Estimation des coefficients finaux par régression logistique multivariée */
%put Calcul des coefficients finaux;

%let model = T1;

%do c = 2 %to &n_comp.;

%let model = &model. T&c.;

%end;

%let merge = Y &model.;

data log_multi;
merge &merge.;
run;

/* Calcul des coefficients c1, c2, etc associé à T1, T2, etc */
proc logistic data = log_multi;
model &Y. (ref = « &ref. ») = &model.;
ods exclude EffectPlot ROCCurve LeverageCPlot LeverageDifChisqPlot InfluencePlots LeveragePlots LeverageDifDevPlot LeveragePhatPlot CBarPlot CPlot DevianceResidualPlot DifChisqPlot DifDeviancePlot LeveragePlot PearsonResidualPlot ModelInfo ResponseProfile Type3 Classification Association ConvergenceStatus FitStatistics GlobalTests OddsRatios ParameterEstimates ClassLevelInfo NObs;
ods output ParameterEstimates = coeffs_finaux;
run;

data Intercept;
set coeffs_finaux;
if Variable = « Intercept »;
rename Estimate = coeffs_finaux;
keep Variable Estimate;
run;

data _null_;
set coeffs_finaux;
if _n_ = 2;
call symput(« coeff »,Estimate);
run;

data coeff_wetoile;
set coeff_wetoile1;
Estimate = Estimate * &coeff.;
Coeffs_finaux = Estimate;
rename Estimate = Coeffs_comp1;
run;

    /* Addition des coefficients Wetoile multiplié par les coefficients C */
%do c = 2 %to &n_comp.;

data _null_;
set coeffs_finaux;
if _n_ = &c.+1;
call symput(« coeff »,Estimate);
run;

data coeff_wetoile&c.;
set coeff_wetoile&c.;
Estimate = Estimate * &coeff.;
rename Estimate = Coeffs_comp&c.;
run;

data coeff_wetoile;
merge coeff_wetoile coeff_wetoile&c.;
Coeffs_finaux = Coeffs_finaux + Coeffs_comp&c.;
run;

%end;

data Coeffs_finaux;
set Intercept coeff_wetoile;
keep Variable Coeffs_finaux;
run;

    /* Sortie des performances */
%put Calcul des prédictions;

proc transpose data = Coeffs_finaux out = Coeffs_appli;
id Variable;
run;

data save;
set save;
cleeee = 1;
run;

data coeffs_appli;
set coeffs_appli;
drop _name_;
cleeee = 1;
%do p = 1 %to &nb_var.;
coeff&p. = %scan(&list_var.,&p.);
drop %scan(&list_var.,&p.);
%end;
run;

 /* Reconstruction de la table pour présenter des résultats propres si aucune variable qualitative */
%if nb_varQ = 0 %then %do;

data &DATA.;
merge save Coeffs_appli;
by cleeee;
drop cleeee;
exp = Intercept;
%do p = 1 %to &nb_var.;
exp = exp + %scan(&list_var.,&p.) * coeff&p.;
%end;
exp = exp(exp);
pred = exp/(1+exp);
if pred >= 0.5 then Ypred = 2;
if pred < 0.5 and pred ne . then Ypred = 1;
keep &list_var. &Y. pred Ypred
run;

%end;

 /* Reconstruction de la table pour présenter des résultats propres si au moins une variable qualitative */
%if nb_varQ > 0 %then %do;

data save;
merge save Coeffs_appli;
by cleeee;
drop cleeee;
exp = Intercept;
%do p = 1 %to &nb_var.;
exp = exp + %scan(&list_var.,&p.) * coeff&p.;
%end;
exp = exp(exp);
pred = exp/(1+exp);
if pred >= 0.5 then Ypred = 2;
if pred < 0.5 and pred ne . then Ypred = 1;
keep pred Ypred;
run;

data &DATA.;
merge &DATA. save;
run;

%end;

/* Sortie des performanes */
proc freq data = &DATA.;
tables &Y. * Ypred;
run;

/* Suppression des tables inutiles */
proc datasets lib = work nolist;
delete class prod save Xt coeff_wetoile coeffs_appli reglin log_multi prodprod prodprod_ajout prodprod_ajout2 biblioQ biblio Y log_univ encours Intercept disj_encours coeffs W biblioEC Xtp Xti calc calcN calcD Wnorm;
%do c = 1 %to &n_comp.;
delete Coeffs_comp&c. Wnorm&c. coeff_wetoile&c. T&c.;
%end;
run;

/* Reset options */
options notes nospool;
ods results;

%mend;

Package et fonction R: https://cran.r-project.org/web/packages/plsRglm/index.html

\bullet Bibliographie:

– PLS-regression: a basic tool fo chemometrics de Svante Wold, Michael Sjöström, Lennart Eriksson

– Data mining et statistique décisionnelle de Séphane Tufféry

– Probabilité, analyse des données et Statistique de Gilbert Saporta

– La régression PLS: théorie et pratique, M. Tenenhaus

– Introduction à la régression PLS de Carole BINARD

– The Elements of Statistical Learning de Trevor Hastie, Robert Tibshirani et Jerome Friedman