<jupyter_start><jupyter_text>Introduksjon I denne notebooken skal vi illustrere hvordan maskinlæring kan brukes på helsejournaldata via et (lite) eksempel.Se tilhørende modul på MittUiB for mer informasjon og motivasjon til problemstillingen: https://mitt.uib.no/courses/21357/pages/lab-1-journaler-og-maskinlaering. Fra data samlet inn fra 100.000 pasienter og tilrettelagt og tilgjengeliggjort av Microsoft https://microsoft.github.io/r-server-hospital-length-of-stay/index.html skal vi forsøke å forutsi liggetiden til pasienter som legges inn på sykehuset. Et slikt system (om det var tilstrekkelig robust og nøyaktig nok) ville vært svært nyttig for ressursfordeling på et sykehus, blant annet for fordeling av personell og for planlegging av utskrivninger ved en avdeling. Ved å undersøke hvilke egenskaper som gir størst forklaringsverdi for systemet kan en også (potensielt) avdekke interessante sammenhenger mellom journaldata og pasienters tilstand.   Setup<jupyter_code>%matplotlib inline
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import sklearn
from pathlib import Path
import seaborn as sns
DATA = Path('../data/hospital')<jupyter_output><empty_output><jupyter_text>Les inn og se på data Datasettet består av tre filer:<jupyter_code>LOS = pd.read_pickle(DATA/'hospital')
DICT = pd.read_excel(DATA/'Data_Dictionary.xlsx')
METADATA = pd.read_csv(DATA/'MetaData_Facilities.csv')<jupyter_output><empty_output><jupyter_text>Vi tar en titt på alle tre:<jupyter_code>LOS.head()
LOS.info()
DICT
METADATA<jupyter_output><empty_output><jupyter_text>Observasjoner hittilDet ser ut som at `LOS` er vår hovedkilde. `DICT` gir oss nyttig informasjon om features i `LOS`, mens `METADATA` har informasjon om de ulike avdelingene data er samlet fra.  > **OBS:** Vi har endret litt på features i vår versjon av `LOS` i forhold til orginalen:  - Vi har slått sammen alle opphold på lengre enn 8 dager,  - Alle `rcount` markert som "5+" er satt til 5.   Plots og videre utforsking Hva er fordelingen av pasientopphold?<jupyter_code>LOS_value_count = LOS['lengthofstay'].value_counts()
LOS_value_count
LOS_value_count.plot.barh(figsize=(12,8))
plt.show()<jupyter_output><empty_output><jupyter_text>Vi ser at datasettet er ganske skjevt fordelt. De aller fleste opphold på sykehuset er korte.  Hvordan er kjønnsfordelingen?<jupyter_code>gender_counts = LOS['gender'].value_counts()
gender_counts.plot.bar()
plt.show()<jupyter_output><empty_output><jupyter_text>Litt flere menn (0) enn kvinner (1).  Hvor lenge oppholder de seg på sykehuset?<jupyter_code># Vi konverterer først "mer enn 8" til tallet 9:
LOS_tmp = LOS.copy()
LOS_tmp['lengthofstay'].loc[LOS_tmp['lengthofstay'] == 'more than 8'] = 9 # Setter alle med "mer enn 8" til 9

LOS_tmp['lengthofstay'] = LOS_tmp['lengthofstay'].astype(int) # Skift til heltall

# Box plot:
ax = sns.boxplot(x='gender', y='lengthofstay', data=LOS_tmp)<jupyter_output><empty_output><jupyter_text>Korrelasjoner Hvordan korrelerer de numeriske features i datasettet til liggetiden? Her var beskrivelsen av features:<jupyter_code>DICT<jupyter_output><empty_output><jupyter_text>Vi undersøker korrelasjonene med 'lengthofstay':<jupyter_code>correlation_matrix = LOS_tmp.corr()
correlation_matrix['lengthofstay'].sort_values(ascending=False)<jupyter_output><empty_output><jupyter_text>Observasjoner - `rcount` er svært korrelert med liggetid. Antall besøk i løpet av siste 180 dager forteller mye om forventet liggetid. - `facid` har også relativt høy korrelasjon. Det kan skyldes at de ulike avdelingene har ulik kapasitet og roller med hensyn på befolkningen de tjener. De fleste har ikke vært på sykehuset tidligere i løpet av siste halvår, men en god del har hyppige besøk:<jupyter_code>LOS['rcount'].value_counts().plot.bar()
plt.show()<jupyter_output><empty_output><jupyter_text>Kapasiteten til de ulike avdelingene varierer:<jupyter_code>METADATA
METADATA.plot(x='Id', y='Capacity', kind='bar')
plt.show()<jupyter_output><empty_output><jupyter_text>Prediksjon Vi velger å angripe problemet med å predikere liggetid som et klassifikasjonsproblem: altså, > Er forventet liggetid 0, 1, 2, 3, 4, 5, 6, 7, 8, eller mer enn 8 dager? Det er også mulig å sette opp en regresjonsmodell for dette problemet. Du må gjerne forsøke det!<jupyter_code>X = LOS.drop('lengthofstay', axis=1)
y = LOS['lengthofstay']<jupyter_output><empty_output><jupyter_text>Hold av et testsett til å evaluere modellen For å simulere den relle situasjonen der vår modell skal brukes til å predikere liggetid i det en ny pasient kommer inn holder vi av et testdatasett for å evaluere modellen. Dette er til for å simulere den reelle situasjonen.<jupyter_code>from sklearn.model_selection import train_test_split<jupyter_output><empty_output><jupyter_text>Vi bruker stratifisert splitting (`stratify=y` nedenfor) for å sikre at fordelingen av liggetider blir den samme i treningsdata som i testdata.<jupyter_code>X_train, X_test, y_train, y_test = train_test_split(X, y, stratify=y, random_state=42)<jupyter_output><empty_output><jupyter_text>Velg modell Vi forsøker en random forest klassifikator:<jupyter_code>from sklearn.ensemble import RandomForestClassifier
rf = RandomForestClassifier(n_estimators=100, random_state=42, n_jobs=-1) 
# n_jobs=-1 for å bruke alle tilgjenglige CPUer på maskinen<jupyter_output><empty_output><jupyter_text>Vi trener modellen:<jupyter_code>rf.fit(X_train, y_train)<jupyter_output><empty_output><jupyter_text>...og beregner accuracy på test-data:<jupyter_code>from sklearn.metrics import accuracy_score
y_pred = rf.predict(X_test)
accuracy_score(y_test, y_pred)<jupyter_output><empty_output><jupyter_text>Vi predikerer altså liggetiden med ca. 67% nøyaktighet.   Evaluering av resultatet Er dette en god nøyaktighet? Hvilke features bruker modellen til sine prediksjoner?  Forvirringsmatrise Som vanlig vil forvirringsmatriser fortelle oss mye om en klassifikators ytelse.<jupyter_code>from sklearn.metrics import confusion_matrix
cm = confusion_matrix(y_test, y_pred)
from utils import plot_confusion_matrix
fig, ax = plt.subplots(figsize=(12,10))
_ = plot_confusion_matrix(cm, classes=np.unique(y_test), ax=ax)<jupyter_output><empty_output><jupyter_text>Feature importance<jupyter_code>importances = rf.feature_importances_
# Find index of those with highest importance, sorted from largest to smallest:
indices = np.argsort(importances)[::-1]
for f in range(X.shape[1]): 
    print(f'{X.columns[indices[f]]}: {np.round(importances[indices[f]],2)}')
f, ax = plt.subplots(figsize=(12,8))
plt.barh(X.columns[indices], np.round(importances[indices],2))
plt.xlabel("Relative importance")
plt.show()<jupyter_output><empty_output><jupyter_text>Vi ser at modellen lener seg veldig mye på `rcount`, altså antall ganger personen har vært på sykehuset siste 180 dager.. Dette er ikke overraskende hvis vi husker tilbake til korrelasjonene vi undersøkte tidligere. Mer om det etterpå. Som vi husker fra Lab 0 bør en ikke stole for mye på feature importances i random forest-modeller. Det er bedre å bruke **permutation importance**  Permutation importance Fra Lab 0 vet vi at ideen bak permutation importance er at om en feature er viktig for prediksjonen så bør tilfeldig omstokking av dataene i tilhørende søyle føre til en drastisk forverring av modellen. Dersom en feature derimot er uviktig vil en slik shuffling ikke bety så mye.La oss underøke:<jupyter_code>import eli5
from eli5.sklearn import PermutationImportance
perm = PermutationImportance(rf, random_state=42).fit(X_test, y_test)
eli5.show_weights(perm, feature_names = X_test.columns.tolist())<jupyter_output><empty_output><jupyter_text>Vi ser også her at modellen bruker `rcount` som sin viktigste feature, med god margin.  > **Din tur!** En hypotese om hvorfor `psychologicaldisordermajor` havner såpass høyt oppe er at personer med dette flagget satt til True er spesielle. Blant annet kan det tenkes at de ofte ligger på avdelingen "Behavioral". Undersøk dette. Undersøk gjerne også andre hypoteser du måtte komme på.  Et naturlig spørsmål: er det slik at om antall besøk øker, så øker også forventet liggetid? Eller motsatt? Eller er det en mer komplisert sammenheng? Dette kan vi (til dels) undersøke ved å bruke våre **partial dependence plots**, som også ble innført i Lab 0.  Partial dependence plots Hva hender med predikert liggetid dersom antall besøk øker?<jupyter_code>from pdpbox import pdp
pdp_goals = pdp.pdp_isolate(model=rf, dataset=X_test, model_features=X_test.columns.tolist(), feature='rcount', n_jobs=-1)
pdp.pdp_plot(pdp_goals, 'rcount', ncols=3)
plt.show()<jupyter_output><empty_output><jupyter_text>Vi ser at for de korte oppholdene øker ikke forventet antall dager med antall sykehusbesøk. For oppholdene på 6 dager eller mer er det derimot en slik sammenheng.  > **Din tur!** Hva slags lidelser (av de som er registrert i datasettet vårt) karakteriserer pasienter som har høy `rcount`?  Feature engineering Vi har observert at visse features i datasettet er viktigere enn andre for vår modell. Noen features beskriver pasientene på enn mer *hensiktsmessig* måte enn andre. Det vil si, med høyere verdi for vår modells prediksjoner av liggetid. > Men det er ingenting som tilsier at dette er de *beste* mulige features for pasientene!  En av de aller viktigste delene av maskinlæring (men også svært undervurdert i mer «teoretiske» vinklinger på feltet) er såkalt **feature engineering**. Maskinlæring handler i bunn og grunn om **å tilnærme funksjoner**. En ønsker å finne en funksjon som er så nært som mulig *fasit-funksjonen*: den som sender mengden av features til en instans til korrekt output. Det vil si, til korrekt klasse om klassifikasjon, til korrekt verdi om regresjon. > Jo mer *komplisert* fasit-funksjonen er jo vanskeligere er det for maskinlæringsmodellen å tilnærme seg denne. Ved å beskrive rådata med nye, bedre egnede features kan en gjøre fasit-funksjonen mindre komplisert. Som et ekstremt eksempel: om vi hadde lagt til innleggingsdato og utskrivningsdato til hver pasient som features i vårt datasett hadde fasitfunksjonen vært veldig enkel: utskrivningsdato minus innleggingsdato er lik antall liggedøgn. Det er stort handlingsrom mellom å gjøre ingenting og å rett og slett legge til fasit som en feature (noe som ikke er nyttig)!Det aller beste er om man kan justere på hvordan datainnsamlingen er gjort. Da kan en legge til helt nye features som en tror kan være nyttige for modelleringen (f.eks. legge til alder som en feature i våre data). Hvis man ikke kan endre på hvilke data som er samlet inn kan man tenke følgende: (i) hente inn eksterne features (f.eks. om vi hadde hatt personnummeret til pasientene kunne vi slått dette opp for å finne mer informasjon), eller (ii) lage nye features ved å kombinere de vi har.  > Ethvert maskinlæringsprosjekt vil ha innslag av disse former for feature engineering! (Deep learning er et slags unntak; mer om det senere).  En ny feature: Antall problemer per pasient Hver pasient er tilordnet en rekke flagg av ulike problemer: ``` 'dialysisrenalendstage', 'asthma', 'irondef', 'pneum', 'substancedependence', 'psychologicaldisordermajor', 'depress', 'psychother', 'fibrosisandother', 'malnutrition', 'hemo' ```<jupyter_code>LOS.head()
DICT<jupyter_output><empty_output><jupyter_text>La oss lage en ny feature som teller antall problemer per pasient. Tanken er at antall problemer en pasient har påvirker liggetiden. Om vi legger dette til som en feature gjør vi det enklere for modellen å plukke opp en slik sammenheng. Vi forsøker:<jupyter_code>issues = ['dialysisrenalendstage',
 'asthma',
 'irondef',
 'pneum',
 'substancedependence',
 'psychologicaldisordermajor',
 'depress',
 'psychother',
 'fibrosisandother',
 'malnutrition',
 'hemo']
# Vi legger til summen av antall problemer som en feature, helt på slutten i vår data frame:
LOS.insert(len(LOS.columns)- 1, 'numberofissues', LOS[issues].astype('int').sum(axis=1))
LOS.head()
LOS['numberofissues'].value_counts()
plt.figure(figsize=(12,8))
LOS['numberofissues'].value_counts().plot.bar()
plt.show()<jupyter_output><empty_output><jupyter_text>> **Din tur!** Er pasienter med mange problemer ofte innlagt? Ligger de lenge?  Lag modell med vår nye feature inkludert<jupyter_code>X = LOS.drop('lengthofstay', axis=1)
y = LOS['lengthofstay']
X_train, X_test, y_train, y_test = train_test_split(X, y, stratify=y, random_state=42)<jupyter_output><empty_output><jupyter_text>Vi trener en ny random forest klassifikator:<jupyter_code>rf = RandomForestClassifier(n_estimators=100, random_state=42, n_jobs=-1) 
rf.fit(X_train, y_train)<jupyter_output><empty_output><jupyter_text>...og beregner accuracy på test-data:<jupyter_code>y_pred = rf.predict(X_test)
accuracy_score(y_test, y_pred)<jupyter_output><empty_output><jupyter_text>Vi predikerer altså liggetiden med ca. 68% nøyaktighet. Littegrann bedre enn forrige modell. La oss se på feature importance og permutation importance for å undersøke hvor god vår nye feature er:<jupyter_code>importances = rf.feature_importances_
# Find index of those with highest importance, sorted from largest to smallest:
indices = np.argsort(importances)[::-1]
for f in range(X.shape[1]): 
    print(f'{X.columns[indices[f]]}: {np.round(importances[indices[f]],2)}')<jupyter_output><empty_output><jupyter_text>Vi ser at modellen bruker `number_of_issues` mer aktivt enn hvert av flaggene.  Hva med permutation importance?<jupyter_code>perm = PermutationImportance(rf, random_state=42).fit(X_test, y_test)
eli5.show_weights(perm, feature_names = X_test.columns.tolist(), top=None)<jupyter_output><empty_output><jupyter_text>Her ser vi virkelig effekten av vår nye feature! Den havnet på en klar andreplass, langt over de andre!  Din tur! Vi har lastet ned et annet klassifikasjonsdatasett og plassert det i katalogen `../data/liver/indian_liver_patient.csv`. Det består av data fra 416 lever-pasienter og 167 ikke lever-pasienter samlet inn i Andhra Pradesh, India. Datasettet er nærmere beskrevet her: https://archive.ics.uci.edu/ml/datasets/ILPD+(Indian+Liver+Patient+Dataset). Fra 10 features er oppgaven å predikere hvorvidt en instans svarer til en pasient eller ikke. Last inn data med Pandas, utforsk datasettet som over, og forsøk å predikere pasientstatus.  **Hint:**<jupyter_code>#DATA = Path('../data/liver')
#df = pd.read_csv(DATA/'indian_liver_patient.csv')<jupyter_output><empty_output><jupyter_text>Ekstramateriale Vår accuracy på 68% er ikke særlig god. Det er mange standard måter å øke ytelse i maskinlæring. La oss nevne noen (dette blir fort teknisk; ikke la deg overvelde av de mange ukjente begreper og ideer):1. Få tak i bedre data. Kanskje har en tenkt feil på problemet under datainnsamlingen? Er data formålsmessig samlet inn?2. Få tak i mer data. "Ekte" data, syntetisert data, eller "lignende" data.3. Lag nye, bedre features via feature engineering. Eventuelt, fjern unyttige features (slike kan ødelegge modellen ved at mer data kreves)4. Preprosesser data bedre. Skalering, transformasjoner.5. Bruk en bedre tilpasset modell. De fleste modeller kan "tunes" ved å velge bedre "hyperparametre". Slik kan de gjøres mer eller mindre kompliserte, noe som er nyttig når en skal tilnærme seg en fasitfunksjon av en gitt kompleksitet. «Everything should be made as simple as possible, but not simpler». Dette henger tett sammen med det som kalles "regularisering", "overfitting" og "underfitting".6. Bruk en bedre egnet modelltype. Kanskje har du valgt en type modell som er for simpel (for eksempel logistisk regresjon på for komplisert datasett). Eller kanskje for komplisert?7. Bruk flere modeller sammen i et "ensemble". Hver modell kan for eksempel stemme på et resultat, og så kan en predikere basert på konsensus. "Wisdom of the crowd". Du kan også la hver modell basere seg på litt ulike features, slik at feilene hver modell gjør er mest mulig uavhengige av hverandre. Det er også mulig å bruke en modell til å lære hvordan en samling av ulike modeller best kan settes sammen. Eller sekvensielt lage modeller som forsøker å predikere feilene gjort av tidligere modeller, for slik å sammen minimere feil.  I vårt tilfelle kan vi gjøre et forsøk på nummer 6: vi kan bruke en annen type modell enn random forest.   XGBoost – En av de kraftigste modellene som finnes! (Sammen med dype nevrale nettverk og «deep learning») Vi skal ikke gå inn på hvordan XGBoost fungerer her, bare bruke modellen som en black box (du kan lese om modellen her eller her eller her om du er spesielt interessert). NB: om du ikke har xgboost installert, kjør følgende celle:<jupyter_code>import sys
!{sys.executable} -m pip install xgboost
from xgboost import XGBClassifier
xgb_clf = XGBClassifier()
xgb_clf.fit(X_train, y_train)
xgb_clf.score(X_test, y_test)<jupyter_output><empty_output><jupyter_text>Umiddelbart er ikke dette særlig imponerende.. Dårligere resultat enn vi hadde med random forest! Men det er typisk for XGBoost! Default-verdiene for hyperparameterene i XGBoost (de som ble listet når vi skrev `.fit` over) er sjelden gode. En må velge parametre som er tilpasset problemet. Dette kan gjøres via såkalt **hyperparametersøk**. Dette er utenfor vårt pensum så vi dropper det her. Jeg har tatt et kort søk etter gode parametre (med `GridSearchCV` om du lurer) og endt opp med følgende. Et større søk gjennom flere alternativer kunne ført til bedre resultater. Se her for tips om hvordan en kan optimere XGBoost (om du er særdeles interessert).<jupyter_code>xgb_clf = XGBClassifier(max_depth=15, min_child_weight=15, random_state=0, n_estimators=300, n_jobs=-1)
xgb_clf
#Advarsel: tar en god del tid å trene denne modellen
xgb_clf.fit(X_train, y_train)
xgb_clf.score(X_test, y_test)<jupyter_output><empty_output><jupyter_text>Modellen har en accuracy på 74.3%. La oss ta en titt på permutation importance for denne modellen. **Advarsel: dette tar veldig lang tid!**<jupyter_code># På grunn av en bug i ELI5 når PermutationImportance brukes sammen med XGBoost 
# må vi først trene modellen på nytt på følgende vis (ikke på en Pandas data frame, 
# men på Numpy array. Dette oppnås med `.values`.):
xgb_clf.fit(X_train.values, y_train.values)
perm = PermutationImportance(xgb_clf, random_state=42)
perm.fit(X_test.values, y_test.values)
eli5.show_weights(perm, feature_names = X_test.columns.tolist())<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.07694403574737071, 'fraction_numerical': 0.007193068497629557, 'mean_word_length': 6.154775828460039, 'pattern_counts': {'":': 0, '<': 125, '<?xml version=': 0, '>': 133, 'https://': 3, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '31', 'dir': 'jupyter-structured-clean-dedup', 'id': '649759', 'license': 'MIT', 'n_tokens_mistral': 6774, 'n_tokens_neox': 6499, 'n_words': 2163, 'path': 'nbs/ELMED219-Lab1-Liggetid.ipynb', 'repo_name': 'arvidl/elmed219-brainage', 'size': '18351'}