[Test statistique 2e année / quasi 1re année] Formation à l'analyse de régression en Python (1) --Qiita est une continuation. Cette fois, nous allons en fait essayer une analyse de régression multiple avec deux ou plusieurs variables indépendantes (variables explicatives) en utilisant Python et envisager une interprétation statistique. Cette fois également, j'utiliserai des bibliothèques telles que numpy et pandas pour écrire moi-même la partie calcul en fonction de la formule.
Pour référence de la partie théorique, je me réfère également à ce qui suit. The University of Tokyo Faculty of Liberal Arts Statistics Department "Introduction to Statistics (Basic Statistics I)" Tokyo University Press (1991) Chapitre 13
J'aurais aimé pouvoir le faire avec les données d'itinéraire JR que j'ai utilisées la dernière fois, mais il m'a semblé difficile de collecter les variables indépendantes, j'ai donc décidé d'utiliser une autre donnée cette fois. Page d'accueil du Bureau des statistiques / Statistiques japonaises 2019-Chapitre 2 Population / Ménage
Les données suivantes pour chaque préfecture sont combinées et traitées au format CSV pour utilisation.
―― 2-2 Population par préfecture et taux d'augmentation / diminution de la population ―― 2-13 Population de jour par préfecture et nombre de personnes travaillant / fréquentant l'école en dehors de leur domicile ―― 2-14 Nombre de personnes entrant et sortant par préfecture ―― 2-16 Nombre de naissances / décès et nombre de mariages / divorces par préfecture
** (07/01/2020) Seules les colonnes de population et de ratio de population étaient décalées d'une ligne, je l'ai donc corrigée. ** **
Préfectures,population,population性比,population増減率,昼夜間population比率,Nombre de sur-migrants,taux de natalité,Mortalité,Taux d'augmentation / diminution naturelle,Taux de mariage,Taux de divorce
Hokkaido,5320,89.1,-0.455,99.945,-6569,6.4,11.8,-5.4,4.5,1.92
Aomori,1278,88.6,-0.965,99.849,-6075,6.3,13.8,-7.5,4.0,1.64
:
:
Okinawa,1626,88.5,0.582,99.968,-1112,11.3,8.4,3.0,5.7,2.44
Comme précédemment, nous l'enregistrerons avec le ** code de caractère en UTF-8 et le nom de fichier ** population.csv.
Ci-dessous, nous avons confirmé l'opération avec Python 3.6.8.
import numpy as np
import pandas as pd
df = pd.read_csv("population.csv")
print(df) #Vérifiez le contenu des données importées
En utilisant ces données, pensons à l'analyse de régression lorsqu'il y a deux variables indépendantes ou plus, c'est-à-dire une «analyse de régression multiple». Il existe différentes combinaisons possibles de variables indépendantes et de variables dépendantes, mais tout d'abord, en termes de confirmation de la méthode, il semble qu'une relation facile à comprendre en ressortira.
Essayons avec.
[^ birthrate]: nombre de naissances et de décès pour 1 000 habitants
Considérons un modèle de régression de $ y = \ beta_0 + \ beta_1 x_1 + \ beta_2 x_2 $, où le taux de natalité est $ x_1 $, le taux de mortalité est $ x_2 $ et le taux de croissance démographique est $ y $. Comme pour une variable, coefficient de régression partielle $ \ beta_0 afin de minimiser la somme des carrés des résidus de chaque point de données (différence entre la valeur observée (valeur réalisée) de la variable dépendante et la valeur théorique (valeur prédite)). , \ beta_1, \ beta_2 $. Chaque point de données (valeur observée) $ (x_1 ^ {(1)}, x_2 ^ {(1)}, y ^ {(1)}), ..., (x_1 ^ {(n)}, x_2 ^ { Pour (n)}, y ^ {(n)}) $, la valeur théorique de la variable dépendante est $ \ hat {y} ^ {(i)} = \ beta_0 + \ beta_1 x_1 ^ {(i)} + \ beta_2 Si x_2 ^ {(i)} $, la somme des carrés des résidus $ L $ est
\begin{align}
L &= \sum_{i=1}^n (y^{(i)} - \hat{y}^{(i)})^2 \\
&= \sum_{i=1}^n (y^{(i)} - \beta_0 - \beta_1 x_1^{(i)} - \beta_2 x_2^{(i)})^2
\end{align}
Ce sera. Pour $ \ beta_0, \ beta_1, \ beta_2 $, supposons que la solution lorsque le coefficient différentiel partiel est 0 est la valeur estimée $ \ hat \ beta_0, \ hat \ beta_1, \ hat \ beta_2 $.
\begin{align}
\hat{\boldsymbol{\beta}} &= (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{y}
\end{align}
($ ^ \ Top $ est supposé représenter la transposition d'une matrice / vecteur). pourtant
\begin{align}
\hat{\boldsymbol{\beta}} &= (\hat\beta_0, \hat\beta_1, \hat\beta_2)^\top \\
\boldsymbol{X} &= \begin{pmatrix}
1 & x_1^{(1)} & x_2^{(1)} \\
\vdots & \vdots & \vdots \\
1 & x_1^{(n)} & x_2^{(n)} \\
\end{pmatrix}
\\
\boldsymbol{y} &= (y^{(1)}, ..., y^{(n)})^\top
\end{align}
Supposons donc que $ \ boldsymbol {X} ^ \ top \ boldsymbol {X} $ ait une matrice inverse. En utilisant ce résultat, essayez de calculer le coefficient de régression partielle
n = len(df)
X = np.c_[np.ones(n), df[["taux de natalité", "Mortalité"]].values]
y = df["Taux d'augmentation / diminution de la population"]
b = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(y)
print(b) # Output: [ 1.00978055 0.10175178 -0.18246347]
Par conséquent, $ \ hat {y} = 1,010 + 0,102 x_1 --0,182 x_2 $ est obtenu comme l'équation de régression de la valeur théorique du taux d'augmentation / diminution de la population de $ \ hat {y} $.
――Si le taux de natalité pour 1 000 habitants augmente de 1, le taux annuel d'augmentation / diminution de la population augmentera de 0,102 point. ――Si le taux de mortalité pour 1 000 habitants augmente de 1, le taux annuel d'augmentation / diminution de la population diminuera de 0,182 point.
Si vous essayez de trouver le coefficient de décision comme dans le cas d'une variable
TSS = ((y - y.mean()) ** 2).sum() #Fluctuation totale
ESS = ((X.dot(b) - y.mean()) ** 2).sum() #Fluctuation régressive
RSS = ((y - X.dot(b)) ** 2).sum() #Fluctuation résiduelle
R_2 = ESS / TSS #Coefficient de décision
print(R_2) # 0.9014307808651625
Il a été demandé qu'il soit d'environ 0,901. Cette valeur est-elle grande ou petite? Je vais l'essayer plus tard si je prends une autre variable indépendante, alors comparez les valeurs.
Les détails sont ignorés dans le livre de référence ci-dessus, mais dérivons la formule des coefficients de régression partielle $ \ beta_0, \ beta_1, \ beta_2 $. Il devrait être possible de l'obtenir en la différenciant partiellement normalement et en la résolvant comme une équation simultanée, mais comme la transformation de formule qui était gênante même avec une variable indépendante devient plus gênante, organisons-la avec une vue claire en utilisant la «différenciation du vecteur».
Résumé officiel "Différenciation par vecteur / Différenciation par matrice" - Qiita
Ici, généralisez le nombre de variables indépendantes à $ p $, et trouvez la formule des coefficients de régression partielle $ \ beta_0, \ beta_1, ..., \ beta_p $.
\begin{align}
\boldsymbol{\beta} &= (\beta_0, \beta_1, ..., \beta_p)^\top \\
\boldsymbol{x}^{(i)} &= (1, x_1^{(i)}, ..., x_p^{(i)})^\top \\
\boldsymbol{X} &= (\boldsymbol{x}^{(1)}, ..., \boldsymbol{x}^{(n)})^\top = \begin{pmatrix}
1 & x_1^{(1)} & \cdots & x_p^{(1)} \\
\vdots & \vdots & & \vdots \\
1 & x_1^{(n)} & \cdots & x_p^{(n)} \\
\end{pmatrix}
\\
\boldsymbol{y} &= (y^{(1)}, ..., y^{(n)})^\top \\
\epsilon^{(i)}
&= y^{(i)} - \hat{y}^{(i)} = y^{(i)} - (\beta_0 + \beta_1 x_1^{(i)} + ... + \beta_p x_p^{(i)} ) \\
&= y^{(i)} - \boldsymbol{\beta}^\top \boldsymbol{x}^{(i)} \\
\boldsymbol{\epsilon} &= (\epsilon^{(1)}, ..., \epsilon^{(n)})^\top \\
&= \boldsymbol{y} - \boldsymbol{X} \boldsymbol{\beta}\\
\end{align}
ça ira. Notez que par souci de simplicité, nous avons ajouté une dimension à $ \ boldsymbol {x} ^ {(i)} $, qui correspond au terme constant. Nous écrirons également le vecteur zéro (un vecteur avec toutes les composantes nulles) sous la forme $ \ boldsymbol {0} $. Définir le coefficient différentiel partiel sur 0 pour $ \ beta_0, \ beta_1, ..., \ beta_p $ signifie $ \ frac {\ partial L} {\ partial \ boldsymbol {\ beta}} = \ boldsymbol {0} C'est la même chose que de mettre $.
Recherchez $ \ boldsymbol {\ beta} $ en utilisant la formule de la page ci-dessus, le cas échéant. Il s'agit d'une version légèrement rédigée de la formule avec le contenu suivant. Introduction aux statistiques-Chapitre 7 | Garaku Taton Chenkan> Trivia Room> Trivia Corner
\begin{align}
L &= \sum_{i=1}^n (\epsilon^{(i)})^2 = \boldsymbol{\epsilon}^\top \boldsymbol{\epsilon} \\
&= (\boldsymbol{y} - \boldsymbol{X} \boldsymbol{\beta})^\top (\boldsymbol{y} - \boldsymbol{X} \boldsymbol{\beta}) \\
&= \boldsymbol{y}^\top \boldsymbol{y} - \boldsymbol{y}^\top (\boldsymbol{X} \boldsymbol{\beta}) - (\boldsymbol{X} \boldsymbol{\beta})^\top \boldsymbol{y} + \boldsymbol{\beta}^\top \boldsymbol{X}^\top \boldsymbol{X} \boldsymbol{\beta} \\
&= \boldsymbol{y}^\top \boldsymbol{y} - 2 \boldsymbol{y}^\top \boldsymbol{X} \boldsymbol{\beta} + \boldsymbol{\beta}^\top \boldsymbol{X}^\top \boldsymbol{X} \boldsymbol{\beta} \\
\frac{\partial L}{\partial \boldsymbol{\beta}} &= \frac{\partial}{\partial \boldsymbol{\beta}} \boldsymbol{y}^\top \boldsymbol{y} - 2 \frac{\partial}{\partial \boldsymbol{\beta}}(\boldsymbol{y}^\top \boldsymbol{X}) \boldsymbol{\beta} + \frac{\partial}{\partial \boldsymbol{\beta}}\boldsymbol{\beta}^\top (\boldsymbol{X}^\top \boldsymbol{X}) \boldsymbol{\beta} \\
&= -2 (\boldsymbol{y}^\top \boldsymbol{X})^\top + (\boldsymbol{X}^\top \boldsymbol{X} + (\boldsymbol{X}^\top \boldsymbol{X})^\top) \boldsymbol{\beta} \\
&= -2 \boldsymbol{X}^\top \boldsymbol{y} + 2\boldsymbol{X}^\top \boldsymbol{X}\boldsymbol{\beta}
\end{align}
Par conséquent, quand il y a une matrice inverse dans $ \ boldsymbol {X} ^ \ top \ boldsymbol {X} $, $ \ frac {\ partial L} {\ partial \ boldsymbol {\ beta}} = \ boldsymbol {0} $ Soit la solution de la valeur estimée $ \ hat {\ boldsymbol {\ beta}} $
\begin{align}
\hat{\boldsymbol{\beta}} &= (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{y}
\end{align}
Il peut être obtenu.
Maintenant que nous avons une équation de régression, considérons statistiquement «la fiabilité de l'équation de régression résultante». Jusqu'à présent, je n'ai parlé que d'algèbre linéaire plutôt que de statistiques.
Jusqu'à présent, l'idée était de trouver une équation de régression optimale à partir des valeurs observées de la variable dépendante. Maintenant, que se passe-t-il si ** considérez la variable dépendante comme une variable stochastique qui suit une certaine distribution de probabilité **? En d'autres termes, même si la valeur de la variable indépendante $ \ boldsymbol {X} $ est la même, si vous faites beaucoup d'observations, vous obtiendrez une valeur d'observation différente $ \ boldsymbol {y} $ à chaque fois. À ce stade, le coefficient de régression partielle $ \ hat {\ boldsymbol {\ beta}} $ obtenu par la formule ci-dessus est également une variable stochastique qui change à chaque observation. ** Nous ne connaissons pas la valeur du vrai coefficient de régression partielle, mais nous pouvons considérer la distribution de probabilité que suit le coefficient de régression partielle. ** L'idée à ce sujet est la même que l'estimation par intervalle de la moyenne de la population de la fourchette de niveau 2 du test statistique (la vraie moyenne de la population ne peut pas être connue et la distribution de probabilité de la moyenne de la population est calculée).
Par exemple, lors de la recherche d'une équation de régression, il peut y avoir des variables indépendantes qui jouent un rôle important dans l'équation de régression et des variables indépendantes qui ont peu d'effet. Si $ \ beta_j = 0 $, la variable indépendante correspondante $ x_j $ n'a aucun effet sur la valeur de la variable dépendante. Par conséquent, lorsque l'équation de régression est obtenue à partir de la valeur observée, l'hypothèse nulle $ \ beta_j = 0 $ et l'hypothèse alternative $ \ beta_j \ neq 0 $ sont définies, et si l'hypothèse nulle est rejetée, la variable indépendante est la régression. Vous pouvez le considérer comme une expression qui explique la variable dépendante. Pour cela, il est nécessaire de considérer la distribution de probabilité de $ \ hat {\ beta_j} $ obtenue à partir des valeurs observées.
Nous supposons ici que l'erreur $ \ epsilon ^ {(i)} $ de chaque donnée suit une distribution normale avec une moyenne de 0, variance $ \ sigma ^ 2 $ indépendamment l'une de l'autre. Puis
On peut dire ça. La moyenne de ce $ \ hat {\ boldsymbol {\ beta}} $ correspond au vrai coefficient de régression partielle $ \ boldsymbol {\ beta} $ (preuve ci-dessous). De plus, $ S_ {ij} $ est défini comme le composant $ (i, j) $ de $ (\ boldsymbol {X} ^ \ top \ boldsymbol {X}) ^ {-1} $ (cependant, $ \ beta_j $ La partie supérieure gauche est le composant (0,0) selon la numérotation), et la distribution de chaque $ \ hat {\ beta_j} $ est $ \ sigma ^ 2 S_ {jj} $ (la preuve sera décrite plus loin). Par conséquent, $ \ hat {\ beta_j} $ a été normalisé à 0 en moyenne et à 1 sur la distribution.
Z_j = \frac{\hat{\beta_j} - \beta_j}{\sqrt{\sigma^2 S_{jj}}}
Suivra la distribution normale standard, mais comme d'habitude, il n'est pas possible de connaître la valeur de la variance d'erreur $ \ sigma ^ 2 $. En fait, il n'y a pas d'autre choix que d'estimer à partir de la valeur observée, mais pour en faire une estimation sans biais, le dénominateur est le nombre de données moins (variable indépendante + 1). Soit le nombre de variables indépendantes $ p $
s^2 = \frac{\sum_{i=1}^n (y^{(i)} - \hat{y}^{(i)})^2}{N-p-1}
Est une estimation non biaisée de $ \ sigma ^ 2 $, qui est la variance d'erreur estimée. Remplacement de $ \ sigma ^ 2 $ par $ s ^ 2 $
t_j = \frac{\hat{\beta_j} - \beta_j}{\sqrt{s^2 S_{jj}}}
Suit une distribution $ t $ avec $ n-p-1 $ degrés de liberté. La raison pour laquelle cela suit la distribution $ t $ est similaire à la "statistique t (10.4)" ci-dessous. [Test statistique 2e année] Résumé des distributions de probabilité qui apparaissent dans les tests hypothétiques pour les populations normales - Qiita
Par conséquent, lorsque vous testez si $ \ beta_j = 0 $ peut être dit,
t_j = \frac{\hat{\beta_j} - 0}{\sqrt{s^2 S_{jj}}} = \frac{\hat{\beta_j}}{\sqrt{(y^{(i)} - \sum_{i=1}^n y^{(i)})^2 S_{jj} / (N-p-1)}}
Calculez la valeur de et trouvez la valeur calculée dans la table de distribution $ t $ avec le degré de liberté $ n-p-1 $ pour obtenir la valeur $ P $. Plus la valeur absolue de la valeur $ t $ est élevée, plus la valeur $ P $ est petite, ce qui signifie que la variable indépendante a un sens.
De plus, la variance d'erreur estimée $ s ^ 2 $ obtenue précédemment représente la partie qui ne peut pas être exprimée par l'équation de régression, on peut donc dire que plus elle est petite (relativement) meilleure est (en gros) l'équation de régression. Faisons le.
L'introduction est devenue longue, mais elle est finalement pratiquée. Avec les données ci-dessus, calculons la statistique $ t $ et la valeur $ P $ pour l'hypothèse nulle $ \ beta_j = 0 $ pour les coefficients de régression partielle correspondant aux deux variables indépendantes. Lors du calcul de la valeur $ P $ dans le test, la table de distribution $ t $ est utilisée, mais ici elle est calculée par le module scipy.stats.
from scipy import stats
p = X.shape[1] - 1 #Nombre de variables indépendantes
S = np.linalg.inv(X.T.dot(X)) # (X^T X)Matrice inverse de
s_2 = ((y - X.dot(b)) ** 2).sum() / (n - p - 1) #Valeur estimée de la dispersion des erreurs
t = b / np.sqrt(s_2 * S.diagonal()) #t statistiques
P = (1 - stats.t.cdf(np.abs(t), n - p - 1)) * 2 #Valeur P
print(s_2) # 0.015396963025934803
print(t) # [ 3.273055 3.94856867 -13.93114362]
print(P) # [0.00207535 0.00028009 0. ]
** Dans le test statistique de niveau 2, il y a toujours une question pour lire le résultat de sortie de R, mais cette sortie est obtenue par ce type de calcul, n'est-ce pas? ** **
En regardant ce résultat, $ \ beta_1 = 0 et \ beta_2 = 0 $ sont rejetés à un niveau de signification de 5%. Et on peut lire que le taux de mortalité fonctionne mieux que le taux de natalité. Par contre, à partir de $ s ^ 2 = 0,0154 $, l'erreur standard du taux d'augmentation / diminution de la population [%] est estimée à $ s = 0,124 $, mais je ne peux pas dire si ce résultat est bon ou non car il n'y a pas de cible de comparaison. ..
Ci-dessous, lisez la valeur attendue du vecteur en tant que vecteur dans lequel les valeurs attendues de chaque composant sont disposées. En supposant que la moyenne de l'erreur $ \ epsilon ^ {(i)} $ de chaque donnée est 0
\begin{align}
{\rm E}[\hat{\boldsymbol{\beta}}] &= {\rm E}[(\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{y}] \\
&= (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top {\rm E}[\boldsymbol{y}] \\
&= (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top {\rm E}[\boldsymbol{X} \boldsymbol{\beta} + \boldsymbol{\epsilon}] \\
&= (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{X} \boldsymbol{\beta} \\
&= \boldsymbol{\beta}
\end{align}
Ainsi, la moyenne des coefficients de régression partielle correspond au vrai coefficient de régression partielle $ \ boldsymbol {\ beta} $ (estimation sans biais). Il existe différentes méthodes pour estimer le coefficient de régression partielle autres que la méthode du carré minimum, mais l'obtention d'une estimation sans biais est un grand avantage de la méthode du carré minimum.
En référence à cette page, j'ai ajouté la notation pour compléter le calcul intermédiaire. Section 2.4 Coefficient de régression partielle et son test | Analyse des données / bases et applications
Par la suite, si vous écrivez $ {\ rm V} [\ hat {\ boldsymbol {\ beta}}] $, cela représentera une matrice distribuée co-distribuée. Autrement dit, le composant $ (i, j) $ est la covariance de $ \ hat {\ beta_i}, \ hat {\ beta_j} $ (le composant $ (i, j) $ est la distribution de $ \ hat {\ beta_i} $ Cependant, c'est une matrice telle que la partie supérieure gauche est la composante (0,0) selon la numérotation de $ \ beta_j $).
\begin{align}
\hat{\boldsymbol{\beta}} &= (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{y} \\
&= (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top (\boldsymbol{\epsilon} + \boldsymbol{X} \boldsymbol{\beta}) \\
&= (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{\epsilon} + (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{X} \boldsymbol{\beta} \\
&= (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{\epsilon} + \boldsymbol{\beta}
\end{align}
Utiliser pour être
\begin{align}
{\rm V}[\hat{\boldsymbol{\beta}}] &= {\rm E}[(\hat{\boldsymbol{\beta}} - {\rm E}[\hat{\boldsymbol{\beta}}])(\hat{\boldsymbol{\beta}} - {\rm E}[\hat{\boldsymbol{\beta}}])^\top] \\
&= {\rm E}[(\hat{\boldsymbol{\beta}} - \boldsymbol{\beta})(\hat{\boldsymbol{\beta}} - \boldsymbol{\beta})^\top] \\
&= {\rm E}[((\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{\epsilon})((\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{\epsilon})^\top] \\
&= {\rm E}[(\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{\epsilon} \boldsymbol{\epsilon}^\top \boldsymbol{X} ((\boldsymbol{X}^\top \boldsymbol{X})^{-1})^\top] \\
&= {\rm E}[(\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top \boldsymbol{\epsilon} \boldsymbol{\epsilon}^\top \boldsymbol{X} (\boldsymbol{X}^\top \boldsymbol{X})^{-1}] \\
&= (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top {\rm E}[\boldsymbol{\epsilon} \boldsymbol{\epsilon}^\top] \boldsymbol{X} (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \\
&= (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \boldsymbol{X}^\top (\sigma^2 \boldsymbol{I}) \boldsymbol{X} (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \\
&= \sigma^2 (\boldsymbol{X}^\top \boldsymbol{X})^{-1} \\
\end{align}
($ \ Boldsymbol {I} $ est une matrice unitaire). Si $ S_ {ij} $ est défini comme le composant $ (i, j) $ de $ (\ boldsymbol {X} ^ \ top \ boldsymbol {X}) ^ {-1} $, chaque $ \ hat { Vous pouvez voir que la distribution de \ beta_j} $ est $ \ sigma ^ 2 S_ {jj} $.
Jusqu'à présent, nous avons utilisé deux variables indépendantes, le taux de natalité et le taux de mortalité, afin d'estimer le taux de changement de population. Ici, nous examinerons si cette variable indépendante est vraiment bonne ou s'il est préférable de combiner d'autres variables indépendantes. Il existe peut-être d'autres facteurs qui affectent le taux de changement de la population en plus du taux de natalité et du taux de mortalité.
Si vous sélectionnez arbitrairement deux variables indépendantes à partir de variables autres que le taux d'augmentation / diminution de la population (et le nom de la préfecture qui est l'échelle nominale), quel type de combinaison est le meilleur?
import itertools
y = df["Taux d'augmentation / diminution de la population"]
vars = [x for x in df.columns if x not in ["Préfectures", "Taux d'augmentation / diminution de la population"]]
for exp_vars in itertools.combinations(vars, 2): #Choisissez deux variables indépendantes
X = np.c_[np.ones(n), df[list(exp_vars)].values]
b = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(y) #Calculer le coefficient de régression partielle
TSS = ((y - y.mean()) ** 2).sum() #Fluctuation totale
ESS = ((X.dot(b) - y.mean()) ** 2).sum() #Fluctuation régressive
R_2 = ESS / TSS #Coefficient de décision
p = X.shape[1] - 1 #Nombre de variables indépendantes
S = np.linalg.inv(X.T.dot(X)) # (X^T X)Matrice inverse de
s_2 = ((y - X.dot(b)) ** 2).sum() / (n - p - 1) #Valeur estimée de la dispersion des erreurs
t = b / np.sqrt(s_2 * S.diagonal()) #t statistiques
P = (1 - stats.t.cdf(np.abs(t), n - p - 1)) * 2 #Valeur P
print("Variable indépendante: ", exp_vars)
print("Coefficient de régression partiel: ", b)
print("Coefficient de décision: ", R_2)
print("t statistiques: ", t)
print("Valeur P: ", P)
print("Valeur estimée de la dispersion des erreurs:", s_2)
print()
Comme il est trop long de sortir tous les résultats, ce qui suit est un extrait de la sortie pour une combinaison avec un coefficient de décision de 0,9 ou plus. ** (07/01/2020) Le résultat a été remplacé en raison de la correction des données. ** **
Variable indépendante: ('population', 'Taux d'augmentation / diminution naturelle')
Coefficient de régression partiel: [1.35396495e-01 2.92466698e-05 1.38354513e-01]
Coefficient de décision: 0.9202813093155153
t statistiques: [ 2.50498414 4.06232118 16.12219163]
Valeur P: [0.01602233 0.00019698 0. ]
Valeur estimée de la dispersion des erreurs: 0.012452424232620031
Variable indépendante: ('Nombre de sur-migrants', 'Taux d'augmentation / diminution naturelle')
Coefficient de régression partiel: [2.36387964e-01 6.30652813e-06 1.43339245e-01]
Coefficient de décision: 0.9235317925311585
t statistiques: [ 6.24974791 4.36740999 18.56543006]
Valeur P: [1.44848286e-07 7.53394644e-05 0.00000000e+00]
Valeur estimée de la dispersion des erreurs: 0.011944683881961054
Variable indépendante: ('taux de natalité', 'Mortalité')
Coefficient de régression partiel: [ 1.00978055 0.10175178 -0.18246347]
Coefficient de décision: 0.9014307808651792
t statistiques: [ 3.273055 3.94856867 -13.93114362]
Valeur P: [0.00207535 0.00028009 0. ]
Valeur estimée de la dispersion des erreurs: 0.015396963025934803
Variable indépendante: ('taux de natalité', 'Taux d'augmentation / diminution naturelle')
Coefficient de régression partiel: [ 1.02915847 -0.08316346 0.18254551]
Coefficient de décision: 0.9030150527461578
t statistiques: [ 3.35453348 -2.39398938 14.07003206]
Valeur P: [0.00164462 0.02099475 0. ]
Valeur estimée de la dispersion des erreurs: 0.01514949250939055
Variable indépendante: ('Mortalité', 'Taux d'augmentation / diminution naturelle')
Coefficient de régression partiel: [ 1.00487249 -0.08024897 0.10178654]
Coefficient de décision: 0.902445538157104
t statistiques: [ 3.30081126 -2.33256484 4.02629601]
Valeur P: [0.00191782 0.02430887 0.0002203 ]
Valeur estimée de la dispersion des erreurs: 0.01523845329397937
Variable indépendante: ('Taux d'augmentation / diminution naturelle', 'Taux de mariage')
Coefficient de régression partiel: [-0.56327258 0.12603311 0.16092148]
Coefficient de décision: 0.9001397075218234
t statistiques: [-1.3456451 7.29460315 2.0734612 ]
Valeur P: [1.85310560e-01 4.23937485e-09 4.40154160e-02]
Valeur estimée de la dispersion des erreurs: 0.015598634589388555
Dans ce cas, est-ce acceptable parce que la combinaison du nombre de sur-migrants et du taux naturel d'augmentation / diminution a un coefficient de détermination élevé et une petite valeur estimée de la dispersion des erreurs? Cependant, le nombre de sur-migrants est réel et non en pourcentage, de sorte que l'impact sur le taux naturel d'augmentation ou de diminution devrait varier en fonction de la population d'origine. Ici, au contraire, ce qui est divisé par la population devrait être appelé le taux de surmigration et être une variable indépendante.
df["Taux de migration excessive"] = df["Nombre de sur-migrants"] / (df["population"] * 1000) * 100
n = len(df)
exp_vars = ("Taux de migration excessive", "Taux d'augmentation / diminution naturelle")
y = df["Taux d'augmentation / diminution de la population"]
X = np.c_[np.ones(n), df[exp_vars].values]
b = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(y)
TSS = ((y - y.mean()) ** 2).sum() #Fluctuation totale
ESS = ((X.dot(b) - y.mean()) ** 2).sum() #Fluctuation régressive
R_2 = ESS / TSS #Coefficient de décision
p = X.shape[1] - 1 #Nombre de variables indépendantes
S = np.linalg.inv(X.T.dot(X)) # (X^T X)Matrice inverse de
s_2 = ((y - X.dot(b)) ** 2).sum() / (n - p - 1) #Valeur estimée de la dispersion des erreurs
t = b / np.sqrt(s_2 * S.diagonal()) #t statistiques
P = (1 - stats.t.cdf(np.abs(t), n - p - 1)) * 2 #Valeur P
print("Variable indépendante: ", exp_vars)
print("Coefficient de régression partiel: ", b)
print("Coefficient de décision: ", R_2)
print("t statistiques: ", t)
print("Valeur P: ", P)
print("Valeur estimée de la dispersion des erreurs:", s_2)
print()
Variable indépendante: ('Taux de migration excessive', 'Taux d'augmentation / diminution naturelle')
Coefficient de régression partiel: [0.20976966 0.75632976 0.10990394]
Coefficient de décision: 0.9580723972439072
t statistiques: [ 7.50186534 8.42827532 14.31096232]
Valeur P: [2.11598050e-09 9.87150361e-11 0.00000000e+00]
Valeur estimée de la dispersion des erreurs: 0.006549283387531298
Par rapport à l'utilisation du nombre de transfert en excès tel quel, le coefficient de détermination et la valeur estimée de la dispersion d'erreur sont meilleurs. À noter que l'on pense que le taux de transfert excédentaire et le taux naturel d'augmentation / diminution ont pour effet d'augmenter la population, de sorte que le coefficient de régression partielle correspondant est certainement une valeur positive.
Qu'en est-il de l'utilisation de trois variables indépendantes?
for exp_vars in itertools.combinations(vars, 3): #Choisissez 3 variables indépendantes
X = np.c_[np.ones(n), df[list(exp_vars)].values]
b = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(y) #Calculer le coefficient de régression partielle
TSS = ((y - y.mean()) ** 2).sum() #Fluctuation totale
ESS = ((X.dot(b) - y.mean()) ** 2).sum() #Fluctuation régressive
R_2 = ESS / TSS #Coefficient de décision
p = X.shape[1] - 1 #Nombre de variables indépendantes
S = np.linalg.inv(X.T.dot(X)) # (X^T X)Matrice inverse de
s_2 = ((y - X.dot(b)) ** 2).sum() / (n - p - 1) #Valeur estimée de la dispersion des erreurs
t = b / np.sqrt(s_2 * S.diagonal()) #t statistiques
P = (1 - stats.t.cdf(np.abs(t), n - p - 1)) * 2 #Valeur P
print("Variable indépendante: ", exp_vars)
print("Coefficient de régression partiel: ", b)
print("Coefficient de décision: ", R_2)
print("t statistiques: ", t)
print("Valeur P: ", P)
print("Valeur estimée de la dispersion des erreurs:", s_2)
print()
Par exemple, vous pouvez trouver une telle combinaison. (Ce n'est pas le meilleur. C'est un exemple qui facilite l'interprétation des résultats.)
Variable indépendante: ('Taux d'augmentation / diminution naturelle', 'Taux de mariage', 'Taux de divorce')
Coefficient de régression partiel: [-0.30054533 0.12984192 0.19858489 -0.25130329]
Coefficient de décision: 0.9109243128862279
t statistiques: [-0.72219064 7.82613187 2.61426129 -2.28169067]
Valeur P: [4.74086687e-01 8.36975600e-10 1.22779282e-02 2.75178580e-02]
Valeur estimée de la dispersion des erreurs: 0.014237611977607154
Plus le taux naturel d'augmentation / diminution et le taux de mariage sont élevés, plus la population est élevée et plus le taux de divorce est élevé, plus la population est basse. C'est un résultat qui semble convaincant. Les trois coefficients de régression partielle (à l'exception du terme constant) sont significatifs de 5%. D'un autre côté, regardons quelques exemples pas si bons.
Variable indépendante: ('taux de natalité', 'Mortalité', 'Taux d'augmentation / diminution naturelle')
Coefficient de régression partiel: [ 1.03883742 -0.2559058 0.17133451 0.35353742]
Coefficient de décision: 0.9034996234109349
t statistiques: [ 3.3482056 -0.68534239 0.4646739 0.96013752]
Valeur P: [0.00169918 0.49680614 0.6445094 0.34235404]
Valeur estimée de la dispersion des erreurs: 0.015424353851135318
Le résultat étrange est que lorsque le taux de natalité augmente, la population diminue, et que le taux de mortalité augmente, la population augmente. En effet, il existe une relation entre les variables indépendantes «taux d'accroissement / diminution naturel = taux de natalité-taux de mortalité», ce qui pose le problème gênant de la colinéarité multiple. Qu'est-ce que la colinéarité multiple? Dictionnaire Weblio
Lorsque la sélection des variables n'est pas effectuée, il vaut mieux ne pas inclure celles qui présentent une forte corrélation entre les variables indépendantes. Si certains d'entre eux ne sont pas indépendants (par exemple, les variables A et B et leur somme C = A + B sont toutes deux incluses), l'analyse échouera. Dans certains cas, le signe du coefficient de corrélation entre chaque variable indépendante et la variable dépendante peut ne pas correspondre au signe du coefficient de régression partielle. (Omis) Cela se produit parce qu'il existe un mélange de variables indépendantes hautement corrélées (lorsqu'une variable annule la partie sur-prédite avec une autre variable). Il y a). Cependant, cela n'est pas pratique lorsque l'on considère la relation causale, il est donc conseillé de rechercher des équations de régression multiples excluant les variables indépendantes dont les signes ne correspondent pas.
Il peut être possible de détecter si une colinéarité multiple se produit ou non dans une certaine mesure, mais il est difficile pour les machines de comprendre l'histoire selon laquelle «plus le taux de natalité est élevé, plus la population n'augmentera pas». Il semble que «les équations de régression multiple excluant les variables indépendantes dont les signes ne correspondent pas» soit en fait assez difficile. L'interprétation du résultat final est là où les humains ont besoin de réfléchir.
La longue analyse de régression est terminée pour le moment. La prochaine fois, nous prévoyons le chapitre 3 «Analyse des composants principaux».
Recommended Posts