Hvorfor $R^2$ ikke er så nyttig

\(R^2\) eller ‘forklaringsgraden’ (engelsk: coefficient of determination) bliver ofte nævnt som et mål for, hvor godt en regressionsmodel passer til data. Der er dog mange problemer med \(R^2\), som gør, at \(R^2\) ikke er så anvendelig som man skulle tro: \(R^2\) siger ikke noget om, hvor godt modellen passer til data, det siger \(R^2\) ikke noget om, hvor god modellen er til at prædiktere, \(R^2\) kan ikke sammenlignes på tværs af datasæt, og det giver heller ikke mening at sammenligne \(R^2\)-værdier, når man transformerer sine data. Ingen af disse problemer bliver løst ved at bruge justeret \(R^2\), og her illustrerer vi, hvorfor \(R^2\) ikke bør bruges, og hvilke andre metoder man i stedet kan bruge for at belyse specifikke egenskaber ved sin model.

I det efterfølgende vil vi bruge nedenstående datasæt om standselængde og hastighed af biler.

Figur 1: Data er fra 1920erne og findes i R som datasættet cars i pakken ggplot2. Den røde linje angiver den bedste rette regressionslinje.

Een mulig model for disse data er en lineær regressionsmodel

\[ y_i = b_1 + b_2 x_i + e_i, \]

hvor \(y_i\) er standselængder, \(x_i\) er hastigheder og \(e_i\) er residualer. Når man fitter denne model med mindste kvadraters metode kan man beregne fittede eller prædikterede værdier

\[ \hat y_i = \hat b_1 + \hat b_2 x_i \]

Hvis \(\bar y\) er gennemsnittet af \(y_i\)erne så er \(R^2\) defineret som

\[ R^2 = \frac{\sum_i(\hat y_i-\bar y)^2}{\sum_i(y_i-\bar y)^2 } \]

Ovenfor er både tæller og nævner positive og tælleren er mindre end nævneren. Derfor er \(0\le R^2\le 1\). Jo tættere de forventede værdier, \(\hat y_i\), er på tæt på de observerede værdier, \(y_i\), desto tættere vil \(R^2\) være på 1.

For data om standselængde får vi nedenstående estimater for parametrene i regressionsmodellen.

Tabel 1: Parameterestimater for regressionsmodellen anvendt på cars datasættet. Som figuren også antyder er der en kraftig statistisk sammenhæng mellem hastighed og standselængde.
term estimate std.error statistic p.value
(Intercept) -17.58 6.758 -2.60 0.012
speed 3.93 0.416 9.46 0.000

Hastigheden er klart signifikant, og for hver stigning i hastighed (i mph) er der en forventet øget standselængde på lige under 4 fod. \(R^2\) (og justeret \(R^2\) som vi kommer tilbage til nedenfor) bliver i disse data

Tabel 2: Tabel 2: \(R^2\) og justeret \(R^2\) for en lineær regressionmodel anvendt på cars datasættet.
r.squared adj.r.squared
0.651 0.644

Der er som sådan ikke noget forkert i begrebet forklaringsgrad eller coefficient of determination. Spørgsmålet er mere, hvor ofte det rent faktisk er forklaringsgraden, som man er interesseret i. Hvis regressionslinjens hældningen nu er tæt på nul, men punkterne ligger tæt omkring denne linje, så vil \(R^2\) være stor, men den reelle sammenhæng er uden nogen praktisk betydning. \(R^2\) angiver, hvor stor en andel af variationen i data som modellen “forklarer”. Problemet ved at betegne \(R^2\) som forklaringsgraden er, hvis ordet “forklaring” bliver opfattet som en kausal sammenhæng. I eksemplet “forklarer” hastigheden 65% af variationen i standselængder, men standselængderne forklarer dermed også 65% af variationen i hastighed. Betegnelsen “forklarer” er derfor uhenigtsmæssig, fordi den antyder, at den ene variabel er årsagen til den anden.

Andre fortolkninger af \(R^2\)

Der findes mere håndgribelige fortolkninger af \(R^2\): Hvis modellen stemmer godt overens med data så skal \(y_i\) og \(\hat y_i\) ligne hinanden, så det vil være naturligt at plotte disse mod hinanden (se nedenfor):

  1. Hvis modellen er god til at beskrive data, så vil man forvente, at de observerede og fittede værdier skal være tæt på hinanden, og korrelationen mellem disse to størrelser dermed vil være tæt på \(1\). Hvis modellen indeholder et konstantled (\(b_1\) i regressionsmodellen ovenfor) så er kvadratet på denne korrelation netop lig med \(R^2\).

  2. Tilsvarende kan man - igen hvis modellen indeholder et konstantled - beregne \(R^2\) som hældningen på regressionslinien, hvor \(\hat y_i\) er respons og \(y_i\) er forklarende variabel.

  3. Endelig kan man i det tilfælde, hvor man ser på en simpel lineær regressionsmodel som ovenfor beregne \(R^2\) som den kvadrerede korrelation mellem \(y_i\)’erne og \(x_i\)’erne.

\(R^2\) siger ikke noget om, hvorvidt modellen er god

Ovenfor har fokus været på, om \(\hat y_i\)’erne ligger tæt på \(y_i\)’erne, men i virkeligheden siger det ikke noget om, hvorvidt den statistiske model er korrekt. Selv hvis man kender den sande model kan man godt have situationer, hvor \(R^2\) er tæt på 0 ligesom at man kan have en forkert model, hvor \(R^2\) er tæt på 1. Hvis vi bruger vores lineære regressionseksempel ovenfor og antager, at vi kender de sande parametre i modellen (\(b_1\) og \(b_2\)) så vi ikke skal estimere dem, så kan vi skrive \(R^2\) som

\[R^2 = \frac{b_2^2 V(X)}{b_2^2 V(X) + \sigma^2},\]

hvor \(\sigma^2\) er variationen omkring regressionslinjen, og hvor \(V(X)\) angiver variationen af \(x\)’erne. Bemærk, at selv hvis vi kender den sande model så kan vi have en værdi for \(R^2\), der vil være tæt på 0, hvis \(\sigma^2\) er stor. Tilsvarende kan vi opnå en større værdi for \(R^2\), hvis vi kan sprede vores \(x\)-værdier så bredt ud som muligt. I vores eksempel med standselængder kan man se, at hvis man splitter data op i intervaller for hastighed, og udregner \(R^2\) i hvert interval så får man betragtelig mindre \(R^2\)-værdier for hvert interval.

Figur 2: \(R^2\) for alle data er 0.65, men når vi udregner \(R^2\) i hvert af de mindre intervaller fra 0-10 mph, 10-20 mph, 20-30 mph så får vi \(R^2\)-værdier som alle er meget mindre fordi forholdet mellem spredningen af punkterne og spredningen af \(x\)-værdierne ændrer sig.

\(R^2\) vokser med stigende modelkompleksitet

Der er et andet notorisk problem med \(R^2\): \(R^2\) vokser jo flere forklarende variable man tager med i modellen uanset om disse variable har en reel sammenhæng eller betydning med udfaldet. Det er ikke så svært at overbevise sig om, at det må være tilfældet, for en model med flere forklarende variable må nødvendigvis være bedre til at beskrive data end en model, hvor nogle af de forklarende variable ikke er med.

Figur 3: I figuren har vi tilpasset to modeller: en ret linje (den røde linje) og et 8.-gradspolynomium (den blå kurve). Den justerede \(R^2\) for modellen med den rette linje er 0.64, mens den justerede \(R^2\) for 8.-gradspolynomiet er 0.65. Til sammenligning er \(R^2\) for 8.-gradspolynomiet 0.70

Hvorfor er det at et 8. grads polynomium giver en absurd model? Vi kan prøve at tænke i retning af prædiktioner af fremtidige standselængder. Den mindste (største) hastighed i datasættet er 4 mph (25 mph). Hvis vi forsøger at lave prædiktioner blot en smule udenfor dette interval går det helt galt (men også indenfor dataområdet bliver prædiktionerne gale, som det fremgår af figuren ovenfor). \(R^2\) bliver som nævnt større, når man tager flere parametre med i modellen, men samtidig risikerer man at ende med modeller, der er absurde. For at reparere på dette bruger man ofte den justerede determinationskoefficient \(\bar R^2\), (adjusted \(R^2\)) i stedet, da denne ikke nødvendigvis stiger, når flere variable tages med.

Estimation af parametre i en lineær model er baseret på at minimere residual summen af kvadrater \[RSS_1 = \sum_i (y_i - \hat \mu_i)^2\]

Lad \(\bar y\) betegne gennemsnittet af alle \(y_i\)’erne og definer

\[RSS_{0} = \sum_i (y_i-\bar y)^2, \quad RSS_{10} = \sum_i (\hat y_i-\bar y)^2\]

Alle disse størrelser er ikke-negative, og det er ikke svært at vise at \[RSS_0=RSS_1+RSS_{10},\] hvilket betyder at \[ R^2=\frac{RSS_{10}}{RSS_{0}} = \frac{RSS_0-RSS_1}{RSS_{0}} = 1-\frac{RSS_{1}}{RSS_{0}} \]

Bemærk, at \(RSS_0\) er et mål for den samlede variabilitet i data (og denne størrelse har intet at gøre med en model). Minimering af $ RSS_1 $ skal derfor være den samme som maksimering af $ RSS_{10} $.

Residualvariansen \(\sigma^2\) estimeres som \(\hat\sigma^2 = RSS_1/(n-p)\) hvis der er \(p\) forklarende variabler i modellen. Under den simple model hvori der kun er det konstante led estimeres residualvariansen som \(\hat\sigma^2=RSS_0/(n-1)\). Disse to forskellige variansestimater indgår i den justerede \(\bar R^2\)

\[ \bar R^2= 1-\frac{RSS_{1}/(n-p)}{RSS_{0}/(n-1)} \]

Bemærk at \(\bar R^2\) vokser ikke nødvendigvis selvom man tilføjer flere parametre. \(\bar R^2\) fortolkes ikke umiddelbart som “den andel af variation af data som modellen forklarer”, og selvom \(\bar R^2\) forsøger at håndtere problemet med at \(R^2\) vokser med antallet af parametre, så løser det stadig ikke problemet med, at vi ikke kan bruge værdien til at sige noget om, hvorvidt modellen er rimelig eller korrekt. Endelig er der intet der forhindrer, at den justerede \(R^2\) bliver negativ.

Figur 4: Plot af \(R^2\) (rød) og justeret \(\bar R^2\) (blå) mod graden af det polynomium, der fittes til data i cars datasættet.

Modelsammenligning baseret på \(R^2\)

Konklusionen ovenfor er, at man skal være meget forsigtig med at sammenligne modeller alene baseret på deres \(R^2\) værdi. Men man kan godt sammeligne nestede modeller baseret på deres \(R^2\) værdier og antallet af parametre i modellen. Lad os fokusere på to modeller: Modellen baseret på et 2. grads polynomium og modellen baseret på et 1. grads polynomium (en ret linie). Lad os kalde dem for \(M_c\) (for complex) og \(M_r\) (for reduceret). \(M_r\) er en delmodel af \(M_c\), hvilket her bare betyder, at man kan tage en parameter i \(M_c\) og sætte til at være \(0\) og så får man \(M_r\). Man siger også at \(M_r\) er nested i \(M_c\). Bemærk, at det kun gælder for nestede modeller, der sammenlignes på det samme datasæt. Der er derfor ikke muligt at bruge \(R^2\) til at sammenligne to modeller, hvor man for eksempel har logaritme-transformeret responset i den ene model, eller hvis man sammenligner \(R^2\) værdier fra to forskellige datasæt.

Hvilken betydning har \(R^2\) så?

Hvis det primært formål med analysen er at beskrive sammenhængen mellem udfaldet og de forklarende variable så er \(R^2\) ikke til megen hjælp. Hvis den statistiske model er korrekt så ændrer \(R^2\) ikke på den måde man skal fortolke resultaterne (parameterestimaterne) på uanset hvilken værdi, som \(R^2\) måtte have. Hvis det primære formål med analysen er at prædiktere udfaldet som funktion af en række forklarende variable, så er der noget information man kan få ud af \(R^2\). I forbindelse med prædiktioner er det nemlig ikke kun interessant at prædiktere en enkelt værdi, men også at angive usikkerheden på prædiktionen. Når \(R^2\) er tæt på 0 så betyder det, at modellen har større usikkerhed, og omvendt, at man laver prædiktioner med lille usikkerhedsmargin, hvis \(R^2\) er tæt på 1. Men så hjælper \(R^2\) heller ikke længere, for vi kan ikke bruge værdien af \(R^2\) til at sige noget om, hvorvidt usikkerheden af prædiktionen er lille nok i en given situation. For eksempel root mean squared error \(\sqrt{\frac1n \sum_{i=1}^n (\hat y_i - y_i)^2}\) Med andre ord vil det alligevel være nødvendigt at udregne endnu et mål for prædiktionsusikkerheden, og hvis man alligevel skal gøre det, så kunne man jo lige så godt have startet med det til at begynde med.