Regresjon med prognoseskΓ₯r for forverring av diabetes.ΒΆ
Denne notatboken ble brukt pΓ₯ Zoom-forelesning nr 2, 28.10.2020, se ogsΓ₯ forelesningsnotater fra timen pΓ₯ https://www.math.ntnu.no/emner/IST100x/ISTx1003/Zoom2InClass20201028.pdf
SpΓΈrsmΓ₯l: hvilke kovariater kan gi forverret prognoseskΓ₯r?
In [1]:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from scipy import stats
import statsmodels.formula.api as smf
import statsmodels.api as sm
import warnings
warnings.simplefilter(action='ignore', category=FutureWarning)# fordi en av plottefunksjonene mΓ₯ oppdateres av utvikler og kommer per i dag med lang FutureWarning som vi ikke vil se pΓ₯
from IPython.core.interactiveshell import InteractiveShell
InteractiveShell.ast_node_interactivity = "last_expr"
Steg 1: Bli kjent med dataene ved Γ₯ se pΓ₯ oppsummeringsmΓ₯l og ulike typer plottΒΆ
Datasettet inneholder flere kovariater enn vi ser pΓ₯ i forelesningen, se i slutten av filen for tilpasning der de ogsΓ₯ er med. Vi ser pΓ₯ age, sex, bmi, ldl og hdl kolesterol.
In [2]:
filurl="https://web.stanford.edu/~hastie/CASI_files/DATA/diabetes.csv" # fra lærebok i statistisk læring
df=pd.read_csv(filurl)
print(df.shape)
df.head() # vi ser vi har 442 observasjoner av 11 varialber - en er repons (prog)
Out[2]:
In [3]:
# boksplott av kategoriske variabler
# kryssplott av kontinuerlige varialber
ax = sns.boxplot(x="sex", y="prog", data=df)
plt.show()
sns.relplot(x='age', y='prog',data = df)
plt.show()
sns.relplot(x='bmi', y='prog',data = df)
plt.show()
sns.relplot(x='ldl', y='prog',data = df)
plt.show()
sns.relplot(x='hdl', y='prog',data = df)
plt.show()
# et fancy plott der alle kontinuerlige er sammen - og der ser vi ogsΓ₯ om kovariatene varierer sammen
sns.pairplot(df, vars = ['prog','age', 'bmi','hdl','ldl'],
diag_kind = 'kde',
plot_kws=dict(alpha=0.4))
plt.show()
Steg 2: Spesifiser en matematisk modellΒΆ
In [4]:
formel='prog~age+sex+bmi+ldl+hdl'
# andre formler kunne vært
#formel='prog~ldl' # for enkel lineær regresjon med bare ldl
#formel='prog~age+sex+bmi+map+tc+ldl+hdl+tch+ltg+glu' # hvis alt i data skulle med
Steg 3: Initialiser og tilpass modellenΒΆ
In [5]:
modell = smf.ols(formel,data=df)
resultat = modell.fit()
Steg 4: Presenter resultater fra den tilpassede modellenΒΆ
Her er det vi skal gjΓΈre alt arbeidet - vi skal forstΓ₯ de ulike delene av utskriften!
In [6]:
print(resultat.summary())
Legg til tullekovariat for Γ₯ se hva som skjer med R^2ΒΆ
In [7]:
np.random.seed(0)
IQ = np.random.normal(100,16,442)
formel2='prog~age+sex+bmi+ldl+hdl+IQ'
modell2 = smf.ols(formel2,data=df)
resultat2 = modell2.fit()
print(resultat2.summary())
print(resultat2.rsquared_adj)
Steg 5: Evaluer om modellen passer til dataeneΒΆ
In [8]:
sns.scatterplot(resultat.fittedvalues, resultat.resid)
plt.ylabel("Residual")
plt.xlabel("Predikert verdi")
plt.show()
# Lage kvantil-kvantil-plott for residualene
sm.qqplot(resultat.resid,line='45',fit=True)
plt.ylabel("Kvantiler i residualene")
plt.xlabel("Kvantiler i normalfordelingen")
plt.show()
Til slutt: tar med alle kovariateneΒΆ
In [9]:
formel='prog~age+sex+bmi+map+tc+ldl+hdl+tch+ltg+glu' # hvis alt i data skulle med
modell = smf.ols(formel,data=df)
resultat = modell.fit()
print(resultat.summary())
sns.scatterplot(resultat.fittedvalues, resultat.resid)
plt.ylabel("Residual")
plt.xlabel("Predikert verdi")
plt.show()
# Lage kvantil-kvantil-plott for residualene
sm.qqplot(resultat.resid,line='45',fit=True)
plt.ylabel("Kvantiler i residualene")
plt.xlabel("Kvantiler i normalfordelingen")
plt.show()
In [ ]: