På baggrund af aktindsigt i korrespondancen mellem DCE, DJ og KU fremlægger vildtbiolog Egon Bennetsen her sin konklusion vedrørende tandsnitsmetodens brugbarhed til aldersbestemmelse af krondyr. Og den er mildt sagt væsensforskellig fra Jægerforbundets!
Tekst og foto: Egon Bennetsen
Samlet konklusion vedrørende tandsnitsmetodens præcisionVedrørende R2:Som nævnt i indledningen, hævder Danmarks Jægerforbund (DJ): ” – Tandsnitmetoden til aldersbestemmelse er det bedste, vi har, og nu har Københavns Universitet også gennemregnet DCE’s resultater, og det kommet frem til de samme resultater.” Henrik Meilby fra Københavns Universitet (KU) siger intet om, at tandsnitsmetoden er det bedste vi har, og da slet ikke, at den er anvendelig.  Henrik Meilby bekræfter, at R2 = 0,92 ved regression på de 37 punkter, hvis man anvender den model, der bruges af statistikprogrammerne ved skæring i 0,0. (Model 2) Men han gør samtidig, uopfordret, opmærksom på, at den derved opnåede R2 – værdi er en helt anden, end hvis man anvender modellen for almindelig lineær regression med fri skæring, som jeg altså mener, er den rigtige at bruge. (Model 1) Og han understreger, at hvis DCE og Jægerforbundet bruger Model 2s høje R2 – værdi, ja så giver det faktisk bedre mening i stedet at bruge tre andre metoder til at vurdere, hvad de 37 punkter viser om tandsnitsmetodens usikkerhed. Vedrørende alder +/- 2 år og alder +/- 10 %:Der er endnu ikke fra DCE fremlagt dokumentation for deres påstande om, at tandsnitsmetoden kan fastlægge krondyrenes alder med en usikkerhed på +/- 2 år henholdsvis +/- 10 %. Kilde 1: DCE-rapport: Bæredygtig kronvildtforvaltning. Videnskabelig rapport nr. 106. 2014 Sunde og Haugaard |
A) Angående tandsnitsmetodens præcision (R2)
Indledning
Danmarks Jægerforbund (DJ) har på deres hjemmeside lagt følgende artikel ud: http://www.jaege’’’rforbundet.dk/om-dj/nyhedsarkiv/2017/fortsaet-med-at-sende-kaeber-ind/:
Her står der bl.a.:
”Derfor har Danmarks Jægerforbund før jul fået Institut for Fødevare- og Ressourceøkonomi ved Københavns Universitet (KU) til at gennemgå DCE’s analyse og datamateriale. Her kommer KU frem til de præcis samme konklusioner som DCE, Aarhus Universitet.
– Tandsnitmetoden til aldersbestemmelse er det bedste, vi har, og nu har Københavns Universitet også gennemregnet DCE’s resultater, og det kommet frem til de samme resultater.”
Det undrede mig, at DJ ikke offentliggjorde korrespondancen med Københavns Universitet (KU) og lod folk selv konkludere, hvad forskeren Henrik Meilby fra KU nåede frem til.
Jeg bad derfor DJ offentliggøre materialet straks, idet jeg i modsat fald ville søge aktindsigt.
Da det ikke skete, har jeg fået aktindsigt i korrespondancen mellem KU, DJ og DCE.
************
DJ er meget specifikke med hensyn til, hvad de ønsker fra KU.
Mads Flinterup (DJ) til Henrik Meilby (Bilag 6):
”Kan jeg bede dig om at køre en simpel lineær regression på de 37 observationer hvor du finder r2 på x=y og altså linjen skærer i 0,0?
I det øvelsen skulle have et omfang der kræve fakturering er jeg klar på det. Det slut produkt jeg skal stå tilbage med er at jeg kan sige at kompetent databehandler ved IFRO KU har gennem regnet data og fået r2 til…”
Og Bilag 7:
”Du er velkommen til at give mig alle de statistiske betragtninger du vil, men jeg skal kun bruger r2 ved skæring i 0,0 i forhold til Jægerforbundets pressemeddelelse i samarbejde med DCE i næste uge.” (Min fremhævelse.)
DJ beder altså om én ting – og understreger, at de kun vil bruge denne ene ting – nemlig: R2 = 0,9189 = 0,92.
Og endelig med skæring i 0,0 ved formlen X = Y, med basis i de 37 punkter. (Der skulle nok have stået:
Y= b X).
Den hopper Meilby ikke på!
Lad mig med det samme slå fast, at Meilby efter min opfattelse har bestræbt sig på – og er lykkedes med – at give en meget udførlig og sober redegørelse for problematikken omkring R2 ved de to modeller med henholdsvis fri skæring med Y -aksen (Y = a + b X) og skæring i 0,0 (Y = b X).
Meilbys gennemgang er samlet i bilag 1
Meilby :
”I den konkrete sammenhæng har R-kvadratet den skavank at man får én værdi hvis man opererer med modellen Y = a + b X og en anden hvis man opererer med modellen Y = b X. Det hænger sammen med at kvadratsummerne som standard ikke beregnes ens i de to tilfælde (totalkvadratsummen beregnes i første tilfælde omkring middelværdien af Y, i det andet omkring 0). Det giver som regel meget forskellige R-kvadrater, så det er muligvis det der er årsagen til jeres diskussion?” (Min fremhævelse).
Nemlig!
I: http://html.smartpub.dk/egaagymnasium/infogeist-eg13-b-ma/matematikb/topic_1-40.html
har jeg fundet denne forklarende figur.
Fig. 1 SAK (sum af kvadrat) størrelser
Afsnit 1 – Om regressionsmodeller
Model 1 – Almindelig lineær regression
Y = a X + b (eller Y = a + b X som Meilby kalder den) er vist her.
a X er linjens hældning og b er linjens skæring med Y -aksen.
Regnet ud på baggrund af de 37 punkter, det drejer sig om (Se bilag 1 i: http://www.netnatur.dk/bombe-i-hjortesagen/ ) får vi med formlen 1 – SAK residual/SAK total en R2 værdi på 0,52.
Hvis vi med samme definition sætter b = 0 og tvinger den røde regressionslinje ned i 0,0 får vi 0,50.
I Bilag 2 Meilby. Regneark. Finder vi denne blå boks:
Det er anskueliggjort på denne figur: (Nogle punkter er skjult af andre. (Se bilag 1 i: http://www.netnatur.dk/bombe-i-hjortesagen/ )
Fig.2
R2 kaldes determinationskoefficienten eller forklaringsgraden. Værdien fortæller hvor stor en del af variationen i Y – værdierne, der kan forklares af variation i X -værdierne.
Havde det f.eks. drejet sig om menneskers højde ud ad x -aksen (den uafhængige variabel) og deres vægt op af y – aksen (den afhængige variabel), ville en værdi på 0,50 eller 50 % betyde, at stigningen i vægt for 50 % vedkommende kunne forklares ved øgning i højde.
De andre 50 % måtte skyldes andre faktorer: kropsbygning, alder, træningstilstand m.v.
I dette tilfælde, er x – værdierne dyrenes kendte alder, og y – værdierne de skønnede værdier på baggrund af tandsnit. Her er der altså ikke ”andre faktorer”, der kan forventes at indvirke på y – værdiernes spredning. Enten passer de skønnede værdier med dyrenes reelle værdier, eller også er der en fejl. Det er altså sort/hvidt. Og det R2 måler er metodens præcision. Ved fuld overensstemmelse mellem kendt og skønnet alder ville alle punkter altså befinde sig på en ret linje, begyndende i 2,2 og sluttende i 14,14.
Inden for naturvidenskaberne er kravet til R2 – beregnet på den viste måde – meget høj.
Normalt er kravet at R2 skal være over 95 %.
Her er den 52 % med fri skæring på y -aksen (Y = 1,27) (eller 50 % hvis vi vælger at beregne linjen med skæring i 0,0).
Konklusionen er altså, at metodens præcision er meget svag med hensyn til at bestemme et givet krondyrs alder.
Bemærk: I figuren giver det visuelt fint mening, at R2 kun ændres marginalt fra 0,52 til 0,50 når man ”tvinger” linjen ned gennem 0,0. Punkterne har jo ikke flyttet sig, og linjen har kun foretaget en beskeden drejning.
Et relevant spørgsmål er derfor, hvilke af de 37 punkter vi skal fjerne for at få en værdi på over 0,95 beregnet ved denne model.
Jeg er af den opfattelse, at alle burde kunne blive enige om at anvende et kriterium, der hedder: Kendt alder +/- 1 år. Altså at man kan acceptere at en 8 årig hjort ved tandsnitsmetoden kan blive vurderet til 7, 8 eller 9 år. Det ville lette vurderingen af metodens gyldighed betydeligt, og det kan enhver forstå!
Hvis vi anvender denne regel på de 37 dyr – og dermed fjerner samtlige de 10 punkter (27 %), der afviger mere end +/- 1 år fra kendt alder, ja får vi netop en værdi af R2 på lige over de krævede 95 % (både med og uden skæring i nul).
Fig. 3 (Husk nogle af punkterne er skjult af andre. (Se bilag 1 i: http://www.netnatur.dk/bombe-i-hjortesagen/ )
Det viser med al ønskelig tydelighed, hvor langt tandsnitsmetoden ligger fra, at kunne anvendes til acceptabel bestemmelse af et krondyrs alder.
Model 2 – Lineær regression med fast skæring i 0,0:
Men hvordan kommer Sunde så frem til en R2 værdi på 0,92?
Det sker ved, at statistikprogrammerne i Excel, (bortset fra en enkelt) og SAS (som Peter Sunde bruger) som standard skifter over til en anden regnemetode, når man ”tvinger” linjen gennem 0,0.
Det forklarer Meilby jo fint:
”I den konkrete sammenhæng har R-kvadratet den skavank at man får én værdi hvis man opererer med modellen Y = a + b X og en anden hvis man opererer med modellen Y = b X. Det hænger sammen med at kvadratsummerne som standard ikke beregnes ens i de to tilfælde (totalkvadratsummen beregnes i første tilfælde omkring middelværdien af Y, i det andet omkring 0). Det giver som regel meget forskellige R-kvadrater, så det er muligvis det der er årsagen til jeres diskussion?” (Mine fremhævelser).
Se Meilby. Regneark. Blå kasse. ovenfor, hvor han viser hvad der sker, når man går fra fri skæring af y – aksen til skæring i 0,0 uden at ændre regnemetode. (Se endvidere bilag 3 og 4 i: http://www.netnatur.dk/bombe-i-hjortesagen/)
Se fig. 1 for forklaring:
Det betyder, at Sunde skifter regnemetode ved at dreje den røde linje ned i 0,0, for så skifter statistik-programmet over til at beregne SAK Total fra punkterne og ned til y = 0 (altså X -aksen) i stedet for fra punkterne til y – værdiernes gennemsnit, som er den grønne linje i figur 1.
I ligningen for R2:
betyder det, (Se Meilby. Regneark. Bilag 2)
at medens SAK Residual kun øges fra 214 til 222, så øges SAK Total fra 444 til ikke mindre end 2733.
Derved bliver brøken meget lille og R2 meget stor, nemlig: 0,92.
Og som Meilby helt korrekt gør opmærksom på:” Det giver som regel meget forskellige R-kvadrater…” Min fremhævelse.
I Bilag 3 skriver Sunde: ”Min erfaring fra lignende analyser jeg foretaget i SAS og andre statistiske softwares (bla. SPSS) er nemlig at R2-værdierne stiger en smule, når man går fra alm. lineær regression til lineær regression med fast skæring for datasæt bestående positive værdier.”
I sandhed en flot bemærkning, og noget af en tilsnigelse, når R2 værdien i tandsnits-metoden stiger fra 0,52 til 0,92 – altså med 77 %.
Og det er så altså forklaringen på, at den beskedne drejning af linjen i fig. 2 ikke blot ændrer R2 – værdi fra 0,52 til 0,50, men i SAS-programmet, som Sunde bruger, ændres fra 0,52 til 0,92!
***************
Afsnit 2 – Forkert af DCE at bruge model 2
Hverken Meilby eller jeg er uenige i matematikken bag Sundes beregning af de 92 %, men som anført i mit notat (http://www.netnatur.dk/bombe-i-hjortesagen/) og uddybet i dette afsnit, så mener jeg, at det er en fejl, når DCE/Sunde bruger model 2, da den giver en ”kunstig” høj R2, der ikke afspejler realiteterne.
I Bilag 4 henleder Sunde opmærksomheden på denne kilde:http://www.quepublishing.com/articles/article.aspx?p=2019170
Jeg skylder ham tak herfor!
Linket indeholder en artikel af Conrad Carlberg, hvorfra jeg i Bilag 5 har vist uddrag.
Conrad Carlberg (www.conradcarlberg.com) is a nationally recognized expert on quantitative analysis and on data analysis and management applications such as Microsoft Excel, SAS and Oracle. He holds a Ph.D. in statistics from the University of Colorado and is a many-time recipient of Microsoft’s Excel MVP designation.
Citater fra Carlberg: (Mine fremhævelser og oversættelser)
1) Om skifte af model, når man sætter konstanten til 0:
“Note: No intercept in model. R-Square is redefined.”
Som Meilby, gør også Carlberg opmærksom på, at der er tale om en helt anden R2 – værdi.
2) Om det at skifte model ved at sætte konstanten til 0 skriver han:
Others, including myself, believe that if setting the constant to zero appears to be a useful and informative option, then linear regression itself is often the wrong model for the data.
Det er altså noget, man skal være varsom med at gøre…
3) Om konsekvenserne for R2 værdien ved at sætte konstanten til 0 og dermed vælge DCEs Model 2:
“I should mention that it’s easy to reach the conclusion that forcing the constant to zero results in a more accurate outcome.”…..
“ This (spurious) (Uægte, falske) conclusion is based on a misunderstanding of what happens mathematically when the constant is forced to zero.”
“Of course, the problem is due to the fact that in omitting the constant, we are redefining what’s meant by the term “sum of squares.” As a result, we’re dismembering (Amputere, sønderlemme) the meaning of the R2.”
Med andre ord – ved at skifte fra Model 1 til Model 2 taber R2 – værdien sin betydning som forklaringsgrad.
***
FORKERT VALG AF MODEL
De to R2 – værdier for regressionen henholdsvis med fri skæring og skæring i 0,0 gav en næsten enslydende værdi på henholdsvis 0,52 og 0,50 beregnet med samme definition (Fig. 2). Derfor har jeg hidtil ikke fundet det umagen værd, at starte en diskussion om det rimelige i at ”tvinge” regressionslinjen igennem 0,0.
I Kilde 1 forklares det, at man har lavet undersøgelsen:
”For at verificere, at metoden også kan anvendes til aldersbestemmelse af krondyr under danske forhold…..”.
Om grunden til at fravælge den almindelige lineære regression til fordel for lineær regression med fast skæring (0,0), bemærkes det:
”Denne sammenhæng blev både modelleret som en ret linjemed en fikseret skæring med y-aksen gennem 0 (svarende til at man altid uden usikkerhed vil kunne bestemme alderen på en kalv).”
Hvad sker der lige her?
Der foreligger en metode, hvor man kan lave en entydig bestemmelse af alderen på kalve og 1 års dyr.
Der sker ved at se på antallet af kindtænder. Har dyret 4 (3 præmolarer og 1 molar) er dyret en kalv (1/2 år) og har dyret 5 (3 præmolarer og 2 molarer), er dyret 1 ½ år.
Det bemærkes, at i tabeller og figurer regnes ikke med ½ år, hvad man nok burde, da det er nærmere middelværdien i jagtsæsonen.
Derved får Sunde et tvivlsomt alibi for at give en kalv værdien 0 og dermed lade linjen starte i 0,0.
Metoder:
- Vi har altså en simpel metode (antal tænder), der skal bruges til nemt at fastslå alderen på ½ og 1 ½ årige dyr.
- Og vi har en anden metode (tandsnit), der tænkes anvendt på krondyr fra 2 (2 ½) år og opefter.
Hele øvelsen går altså ud på, at give en fordomsfri (unbiased) vurdering af, hvorledes den skønnede alder ved tandsnit af 37 krondyr i aldersintervallet 2 – 14 år stemmer overens med den faktiske alder.
En sådan fordomsfri vurdering får man jo netop ikke, når man forlods låser linjens ene endepunkt fast i 0,0.
Ovenikøbet med en værdi fra en anden metode og et andet alderssegment.
Specielt når det som her drejer sig om metodetest, er det vigtigt ikke at lade ”kendt viden”- uanset rigtigheden heraf – influere på resultatet.
Dette vil kunne sløre for eventuelle strukturelle fejl ved metoden!
Dropper man denne fastlåsning, får man den linje, der bedst passer til de 37 punkter evalueringen omhandler.
Og det giver ligningen Y = 0,84 + 1,27 med en R2 – værdi på 0,52.
Med andre ord, en linje med hældningskoefficient på 0,84 og skæring med Y -aksen i 1,27.
Samt en chokerende lav R2 – værdi på 0,52!
Og vel at mærke en værdi, der forklarer hvor stor en del af variationen på y -værdierne, der kan forklares af variation i X – værdierne.
I det konkrete tilfælde udtrykt som metodens præcision.
(Som vist i Meilbys Regneark (Bilag 2) er alle modeller enige om denne metode og denne værdi for simpel lineær regression med fri skæring af Y – aksen).
Konklusion på valg af model: (Se fig. 2).
Hvis vi vælger Model 1 med simpel lineær regression uden skæring i nul får vi en R2 – værdi på 0,52 = 52 %.
Det er altså determinationskoefficienten (forklaringsgraden) for den viste regression.
Deri er alle modeller, Meilby – og forhåbentlig også Sunde – enige.
Hvis vi skifter over til Model 2 ved at ”tvinge” linjen i 0,0, får vi en R2 – værdi på 0,92, som altså er ”kunstig” høj og dermed ikke kan anvendes som forklaringsgrad på tandsnitsanalysens præcision.
****
Afsnit 3
3 alternative metoder til vurdering af tandsnitsmetodens præcision
Meilby synes fuldt bevidst om svagheden ved den høje R2 – værdi på 0,92 i Model 2.
I Bilag 1 skriver han:
”Ja, helt perfekt er der ikke noget der er. Mål for usikkerhed findes desværre også i mange former. Frem for R-kvadratet ville jeg nok fokusere mere på spredningen omkring regressionen (residualspredningen / restspredningen) eller et prædiktionsinterval omkring regressionen, for deri ligger der jo et mere reelt udtryk for hvor meget et enkelt dyrs alder typisk afviger fra det der forventes på baggrund af regressionen. Et andet alternativ er helt at undlade at bruge regression og simpelthen kigge på skønnenes afvigelser fra de sande aldre – og deres fordeling. Det kan muligvis være lige så informativt.” (Min fremhævelse).
Mejlby fortrækker altså – frem for R2 værdien på 0,92 – at se på:
- Meilbys første pointe:
(Se fig. 1). Fokusér på den lodrette afstand fra de enkelte punkter til regressionslinjen (rød). (Residualspredningen).
Fig. 4 Fra Sundes SAS beregning.
Residualer: Den lodrette afstand fra de enkelte punkter til regressionslinjen. (Her 0,0 linjen). Bemærk: (Nogle punkter er skjult af andre. (Se bilag 1 i: http://www.netnatur.dk/bombe-i-hjortesagen/ )).Her ses tydeligt, at selvom gennemsnittet ligger tæt på regressionslinjen, så dækker det over en voldsom spredning på begge sider af linjen.
2) Meilbys anden pointe:
Fokusér på prædiktionsinterval omkring regressionslinjen.
Fig. 5 Fra Sundes SAS beregning.
Prædiktion handler om at ”forudsige” nye observationer. 95%-prædiktionsintervallet viser det interval, hvor 95 % af eventuelle nye observationer kan forventes at placere sig.
Med tandsnitsmetodens svage præcision kræver det altså et prædiktionsinterval på mere end +/- 5 år, for at vi kan forvente en 95 % sandsynlighed for, at et nyt punkt vil ligge i intervallet. (De punkterede linjer).
Man kan altså kun ”forudsige”, at stiller vi med et dyr, som vi ved er 8 år, ja, så er der 95 % sandsynlighed for, at den ved tandsnit vil blive vurderet til en alder mellem 3 og 13 år!
Det viser metodens rystende ringe sikkerhed.
3) Meilbys tredje pointe:
”Et andet alternativ er helt at undlade at bruge regression og simpelthen kigge på skønnenes afvigelser fra de sande aldre – og deres fordeling. Det kan muligvis være lige så informativt.”
Nemlig!
Alene ved at se på punkternes spredning i figuren, burde DJ og DCE for længst have droppet enhver forestilling om at bruge tandsnit til aldersvurdering af danske krondyr!
Især hvis de – som forventet at Vildtforvaltningsrådet – vil bruge metoden til at udskille små alderssegmenter i formodet begrænsede geografiske områder. Der er således intet i DCEs beregninger, der viser, at tandsnits-metoden med nogen acceptabel grad af sikkerhed kan fastlægge, hvorvidt 5 % af kæberne fra f.eks. Thy kommer fra hjorte, der er 8 år eller ældre.
************
Samlet konklusion vedrørende R2:
Som nævnt i indledningen, hævder DJ:
”- Tandsnitmetoden til aldersbestemmelse er det bedste, vi har, og nu har Københavns Universitet også gennemregnet DCE’s resultater, og det kommet frem til de samme resultater.”
Henrik Meilby fra KU siger intet om, at tandsnitsmetoden er det bedste vi har, og da slet ikke, at den er anvendelig.
Henrik Meilby bekræfter, at R2 = 0,92 ved regression på de 37 punkter, hvis man anvender den model, der bruges af statistikprogrammerne ved skæring i 0,0. (Model 2)
Men han gør samtidig, uopfordret, opmærksom på, at den derved opnåede R2 – værdi er en helt anden, end hvis man anvender modellen for almindelig lineær regression med fri skæring, som jeg altså mener, er den rigtige at bruge. (Model 1)
Og han understreger, at hvis DCE og Jægerforbundet bruger Model 2s høje R2 – værdi, ja så giver det faktisk bedre mening i stedet at bruge tre andre metoder til at vurdere, hvad de 37 punkter viser om tandsnitsmetodens usikkerhed.
B) Angående DCEs vurdering af tandsnitsmetoden (+/- 2 år)
Se punkt 2 i: http://www.netnatur.dk/bombe-i-hjortesagen/
Her hævdes det (Kilde 1) at:
”I et referencemateriale på 37 individer med kendt alder (mærket som kalve eller 1-årige), kunne de enkelte dyrs alder uden systematiske fejlkilder og med en høj grad af forklaret variation (R2 = 97 %) estimeres vha. antallet af vækstlinjer i tandmateriale, om end med en usikkerhed på ± 2år for det enkelte individ.” Min fremhævelse.
Sunde gør altså i rapporten gældende, at tandsnitsmetoden kan fastlægge den faktiske alder med en usikkerhed på +/- 2 år.
Om dette hævder Sunde nu i Bilag 3: ”Beregner man en gennemsnit og SD på denne fordeling, er gennemsnittet 0,03 (svarer til at hældningen på regressionslinjen var 0,97) og SD = 2,45. …. Det var dette som i rapporten vidst blev løseligt formuleret som ”inden for +/- 2 år” ).”
Nu mener han pludseligt, at det der ”vidst blev løseligt formuleret som ”inden for +/- 2 år” i virkeligheden betyder standardafvigelsen på en fordeling, der viser størrelsen på afvigelserne mellem opgivet alder og alder estimeret ud fra tandsnit!!
Fig. 6 (Sunde Bilag 3)
Det må siges, at være en ganske kreativ tolkning. Både fig. 4 og fig. 6 viser tydeligt, at 7 punkter ligger mere end +/- 2 år fra den sande alder.
Udsagnet om de +/- 2 år er derfor en grov fejl, idet DCE derved får tandsnitsmetoden til at fremstå med en nøjagtighed, der ikke er belæg for.
C) Angående DCEs vurdering af tandsnitsmetoden (+/- 10 %)
I: Modeller for måling af udviklingen i andelen af ældre hjorte i danske krondyrbestande. Notat fra DCE – Nationalt Center for Miljø og Energi Dato: 13. maj 2016, Peter Sunde
hævdes det:
”Aldersbestemmelse vha. tandsnit må betragtes som den mest præcise og eksakte aldersbestemmelsesmetode (+/- 10 % af reel alder) (Sunde & Haugaard 2014))”
Der er ikke i materialet, jeg har fået aktindsigt i nogen bemærkninger til dette.
Der gælder fortsat at 14 af de 37 punkter ikke lever op til påstanden.
Dette er en grov fejl, idet DCE derved får tandsnitsmetoden til at fremstå med en nøjagtighed, der slet ikke er belæg for.
For nærmere om dette, se punktet i http://www.netnatur.dk/bombe-i-hjortesagen/
*****************
Bilag 1:
Fra: Henrik Meilby <heme@ifro.ku.dk>
Sendt: 16. december 2016 18:10:14
Til: Mads Flinterup
Cc: Peter Sunde; Niels Søndergaard
Emne: RE: S1921 har behov for hjælp!
Hej Mads
Ja, helt perfekt er der ikke noget der er. Mål for usikkerhed findes desværre også i mange former. Frem for R-kvadratet ville jeg nok fokusere mere på spredningen omkring regressionen (residualspredningen / restspredningen) eller et prædiktionsinterval omkring regressionen, for deri ligger der jo et mere reelt udtryk for hvor meget et enkelt dyrs alder typisk afviger fra det der forventes på baggrund af regressionen. Et andet alternativ er helt at undlade at bruge regression og simpelthen kigge på skønnenes afvigelser fra de sande aldre – og deres fordeling. Det kan muligvis være lige så informativt.
I den konkrete sammenhæng har R-kvadratet den skavank at man får én værdi hvis man opererer med modellen Y = a + b X og en anden hvis man opererer med modellen Y = b X. Det hænger sammen med at kvadratsummerne som standard ikke beregnes ens i de to tilfælde (totalkvadratsummen beregnes i første tilfælde omkring middelværdien af Y, i det andet omkring 0). Det giver som regel meget forskellige R-kvadrater, så det er muligvis det der er årsagen til jeres diskussion?
Jeg skal gerne regne lidt på jeres data, men inden jeg rører dem skal jeg lige bede dig/jer kigge en ekstra gang på to observationer der ser meget mærkelige ud. Det drejer sig om ID=2010A13 der er skønnet til 14 år men tilsyneladende har en alder på 5 år og ID=2010A14 der er skønnet til at være 15 år men tilsyneladende kun er 7. Er du sikker på at der ikke er tale om tastefejl eller lignende? De falder helt uden for de øvrige afvigelsers fordeling. Eller er det fx noget med at der kan dannes ekstra ringe i forbindelse med graviditet og diegivning? (begge er af hunkøn).
Dbh Henrik
Fra: Henrik Meilby [mailto:heme@ifro.ku.dk]
Sendt: 20. december 2016 16:52
Til: Peter Sunde <psu@bios.au.dk>; Mads Flinterup <mf@jaegerne.dk>
Cc: Niels Søndergaard <nis@jaegerne.dk>
Emne: RE: S1921 har behov for hjælp!
Hej Mads, Peter m.fl.
Godt, tak for de uddybende kommentarer. Så ved jeg hvad problemet er.
Mads, lad mig starte med definitionen af det almindelige R-kvadrat (blot for en ordens skyld). Det er jo et relativt mål for tilpasning som er sådan indrettet at hvis modellen ikke ”forklarer” noget (i gåseøjne fordi det med forklaringskraften ikke nødvendigvis har noget med en egentlig forklaring at gøre), så er R-kvadratet 0, mens det i tilfælde af at modellen så at sige ”forklarer” alt (passer 100% til observationerne) har værdien 1. At modellen ikke ”forklarer” noget falder for lineær regression (Y = a + b*X) sammen med at hældningsparameteren er 0. Det svarer altså til, at den afhængige variabel kan beskrives ved sit gennemsnit (Y = a). Ved bestemmelse af R-kvadratet benyttes den kvadratiske variation omkring gennemsnittet af den afhængige variabel derfor som reference, dvs. som nævner i den brøk af kvadratsummer som R-kvadratet er beregnet ved hjælp af: Hvis SAK(Y) = SUM((Y-Ygennemsnit)^2) og SAK(residual) = SUM((Y-Yprædikteret)^2), beregnet over alle observationerne, så er R-kvadratet = (SAK(Y) – SAK(residual))/SAK(Y) = 1 – SAK(residual)/SAK(Y). Men hvad nu hvis der er tale om en model hvor skæringen forlods har fået værdien 0. Denne model er udtrykt ved Y = b * X, og hvis denne model ikke forklarer noget og altså har hældningsparameteren b = 0, så ender vi med en model der hedder Y = 0. Som reference ved beregning af R-kvadratet benyttes derfor i dette tilfælde som standard den kvadratiske variation omkring 0. Godt så… Når man beregner en værdis afvigelse fra 0 får man værdien selv, så det der i dette tilfælde kommer til at stå i nævneren i R-kvadratet er simpelthen summen af observationernes kvadrater. Det gør at R-kvadraterne for de to modeller bliver ret forskellige.
Godt. Så er der Excel… I er ikke de første der har haft bøvl med Excel’s statistiske funktioner. Som Peter siger er Excel ikke just statistisk software, så der er meget man ikke kan se nogen steder (og i nogle versioner er der også fejl her og der). Netop det problem som I er stødt ind i har jeg jævnligt studerende der bøvler med. Sagen er at Excel har mindst fire forskellige standardredskaber til bestemmelse af en lineær regression. Til bestemmelse af parametrene a og b i modellen Y = a + b * X kan man benytte funktionerne [1] SLOPE()/STIGNING() og INTERCEPT()/SKÆRING(). Til bestemmelse af modeller på formen Y = a + b1*X1 + … + bn*Xn (altså lige fra en helt simpel til en meget stor model) kan man benytte funktionen [2] LINEST()/LINREGR() eller [3] regressionsanalyse-redskabet i data-analyse add-in pakken (som følger med de fleste installationer men skal aktiveres før den er tilgængelig). Endelig kan man til bestemmelse af parametrene i den simple funktion Y = a + b*X benytte [4] trendlinje-redskabet i Excel’s X-Y-plots. Kun de tre sidste metoder beregner R-kvadratet, så i det vedhæftede regneark har jeg vist hvad der sker når man benytter hver af disse til bestemmelse af parametre og R-kvadrater med udgangspunkt i krondyr-datasættet med og uden de to skæve observationer. For overblikkets skyld har jeg indsat analyser både for den almindelige model Y = a + b*X (til venstre) og for den model som det hele drejer sig om, Y = b*X (til højre). Som I kan se til højre er parameterestimater og R-kvadrater i overensstemmelse med Peters resultater beregnet af SAS når man benytter metoderne [2] og [3]. Men hvis man bruger metode [4], trendlinje-redskabet i X-Y plottene, så får man ganske vist det samme hældningsestimat (0.976 med de to skæverter og 0.936 uden), men R-kvadraterne bliver helt anderledes end ved brug af metoderne [2] og [3]. Årsagen kan I se lidt til venstre hvor jeg har indsat en blå boks. Her har jeg beregnet R-kvadratet for den almindelige model (Y = a + b*X) og modellen uden skæring (Y = b*X) ved hjælp af den definition som man normalt anvender for den almindelige model, dvs. med den kvadratiske variation omkring middelværdien af Y som reference. Som det fremgår bliver resultaterne 0.518 for den almindelige model, hvilket alle metoderne [2]-[4] er enige om, og 0.501 for modellen Y = b*X, hvilket lige præcis er hvad trendlinjeredskabet leverer. Årsagen til miseren er altså at graf-redskabet bruger samme definition af R-kvadratet, uanset om der er tale om den ene eller den anden model.
Så kan man naturligvis tænke længe over hvorfor… Det kan selvfølgelig være en teknisk smutter, en forglemmelse, en tanketorsk eller en regulær fejl. Trendlinjeredskabet tillader at man sætter skæringen til præcis det man vil (ikke nødvendigvis 0), så man har også mulighed for at bestemme hældningen i en model der fx hedder Y = 12 + b*X ó (Y-12) = b*X. Det burde have ført til at man havde ændret definitionen af R-kvadratets nævner til kvadratsummen beregnet omkring den fastsatte værdi af skæringen, men det har man altså ikke gjort. Hvorfor? Tja, som Peter siger skal man nok snakke med en af Microsofts programmører (og nok en der efterhånden er lidt oppe i årene men har en god hukommelse, for det redskab har været der ganske længe).
Peter, du skriver at data er tilgængelige for alle interesserede. Er det ok med dig/jer hvis jeg på et tidspunkt bruger dem i en øvelse for mine studerende? De er ganske illustrative.
Jeg håber at det løste problemet. Glædelig jul og godt nytår!
Dbh Henrik
Bilag 2: Meilby Regneark
Bilag 3: Mail fra Sunde til Søndergaard (DJ):
Fra: Peter Sunde [mailto:psu@bios.au.dk]
Sendt: 13. december 2016 17:32
Til: Niels Søndergaard <nis@jaegerne.dk>
Cc: Mads Flinterup <mf@jaegerne.dk>; Claus Lind Christensen <clc@jaegerne.dk>; Lars Haugaard <laha@bios.au.dk>; Aksel Bo Madsen <abm@bios.au.dk>
Emne: Data og resultater for analyse af sammenhæng mellem opgivet og estimeret alder ud fra tandssnit
Kære alle,
Her lidt kommentarer til de analyser og data som ligger til grund for den så berømte Fig.5 i Sunde&Haugaard 2014. Disse data og resultater er tilgængelige for alle som måtte være interesseret da jeg/DCE hverken har ønske om eller interesse i at skjule hvad der ligger på de estimater vi kommer frem til. I kan derfor frit videredistribuere denne email til jeres medlemmer.
For god ordens skyld vil jeg dog sige at jeg forbeholder mig retten til ikke at besvare alle mulige opfølgende spørgsmål af den simple grund at min arbejdstid er knap og må prioriteres mellem mange andre vigtige opgaver end ufinansieret rådgivning om krondyr. Jeg vil understrege, at der i dette ikke ligger noget ønske om ikke at ville forklare data eller resultater til bestemte personer, blot at jeg ikke nødvendigvis har tid til at besvare alle mulige henvendelser altid.
Rådata (som også blev sendt til Egon) findes i excelregnearket.
I filen ”RTF – linear regression.rtf” finder I det oprindelige output fra SAS fra min første analyse af alle 37 dyr med angivet alder.
I filen ”Tandsnit_kendt_alder_uden 2 outliers.rtf” finder i outputtet af den samme analyse (foretaget her i eftermiddags), hvor jeg har taget de to ”skæverter” (outliers) ud som bestod af dyr som var skudt henholdsvis 8 og 9 år over den opgivne alder.
Kommentarer til resultater af den første analyse (n=37 dyr) som blev publiceret i Sunde&Haugaard 2014 (RTF – linear regression.rtf):
R2-værdien på 0.92 som kan læses i outputtet, er defineret som SS(model)/(SS(model) + SS(error). Disse tal fremgår forneden (også taget fra output fra SAS).
At R2-værdien er anderledes når den beregnes i excel er jeg først blevet gjort bekendt med fornylig. Hvorledes R2 beregnes i Excel kan jeg af gode grunde ikke forklare, al den stund der ikke gives yderligere informationer af programmet. Her må man nok kontakte Microsoft for svar (.. jeg har desværre ikke email-adressen på Bill Gates, men det kan være at der andre i Microsoft Inc. som måske kan hjælpe på forespørgsel). Hvis vi skal GÆTTE på hvad forskellen består i, kan jeg kun komme på den forklaring, at excel beregner noget andet og måske har en formelfejl. Den sidste formodning begrunder jeg i, at R2-værdien for simpel regression (y = a + bx), hvilket er den gængse regressionsformel, er den samme for excel og SAS. Min erfaring fra lignende analyser jeg foretaget i SAS og andre statistiske softwares (bla. SPSS) er nemlig at R2-værdierne stiger en smule, når man går fra alm. lineær regression til lineær regression med fast skæring for datasæt bestående positive værdier. Det faktum at R2 i excel falder ved denne ændring af modellen indikerer for mig at den indbyggede formel for R2 (som beregnes forskellig for de to typer af regressionsmodeller) måske ikke ændres i excel. Men dette er dog rene gisninger så længe man ikke kan se hvordan excels (som jo ikke er noget statistiksoftware) beregninger er foretaget. I sådanne tilfælde af forskellige resultater, vil jeg til hver en tid vælge at stole på SAS frem for på excel.
Hvad der måske er nok så interessant, er størrelsen på afvigelserne i mellem opgivet alder og alder estimeret ud fra tandsnit.
Disse fordeler sig som følger:
Beregner man en gennemsnit og SD på denne fordeling, er gennemsnittet 0,03 (svarer til at hældningen på regressionslinjen var 0,97) og SD = 2,45. …. Det var dette som i rapporten vidst blev løseligt formuleret som ”inden for +/- 2 år” ). At der er observationer som afviger mere end dette ses imidlertid også, herunder to uskønne skæverter med afvigelser på 8 og 9 år. Lars og jeg kan ikke give nogen god forklaring på hvorfor disse to dyr ligger SÅ langt fra den angivne værdi, da så store afvigelser ikke er noget vi normalt ser indikeret i forhold til fx tandslid (kan vi sige nu efter en del flere års erfaring med praktisk håndtering af kæber).
Hvis man tænker den kætterske tanke at de to dyr (som stammede fra Oksbøl i 1970erne [tror jeg det var] og var opbevaret i Zoologisk Museums magasiner) kan være have haft forkert angivet alder (ombytninger eller fejlmærkning kan naturligvis være sket, selv på Oksbøl statsskovdistrikt i 1970erne) og derfor tager dem ud af analysen, opnås en gennemsnitlig afvigelse mellem opgivet og estimeret alder på 0,5 år og en SD omkring dette på +/- 1,4 år. Kører man samme regressionsanalyse som før, men blot uden de to meget afvigende observationer, opnås en R2 på 97% (omtrent det samme som for analysen i Fig. 6) og en hældning på 0,94 med 95% konfidensinterval fra 0,88 til 0,99, dvs. lige akkurat statistisk signifikant lavere end 1.
Hvad man vælger at tro på er rigtigt og godt nok, er i denne forbindelse er i sidste ende en trossag. Resultaterne er som resultaterne er. Og disse og analyserne er hermed fremlagt.
Ind til vi evt. er i stand til at foretage en ny kalibrering mellem estimeret alder ved tandsnit og kendt alder fx på dyr fra Thy, er dette det bedste vi har. Det er min faste faglige vurdering at den nu anvendte metode, om end ikke fejlfri, er god nok til at kunne estimere bestandes omtrentlige alderssammensætning med tilstrækkelig sikkerhed til at man kan foretage pragmatiske vurderinger af om målsætninger for hjortes gennemsnitlige, årlige overlevelse er omtrentligt opfyldt. Det er i den forbindelse fuldt muligt (og akademisk meget interessant) at inkorporere effekterne af den her estimerede usikkerhed omkring aldersbestemmelser af enkeltindivider på de samlede aldersfordelinger. Det vil dog alt andet lige betyde at de aldersspecifikke overlevelsesrater estimeret ud fra de tilsyneladende alder-ved-død-fordelinger vil blive beregnet til at blive LAVERE end tilfældet er nu. Det gør vi imidlertid gerne fremadrettet al den stund vi vel alle er interesseret i de mest muligt retvisende estimater for de forvaltede bestandes alderssammensætning og overlevelse J
Med venlig hilsen
Peter
Bilag 4:
Fra: Peter Sunde [mailto:psu@bios.au.dk]
Sendt: 14. december 2016 19:43
Til: Niels Søndergaard <nis@jaegerne.dk>; Mads Flinterup <mf@jaegerne.dk>; Claus Lind Christensen <clc@jaegerne.dk>
Cc: Aksel Bo Madsen <abm@bios.au.dk>; Lars Haugaard <laha@bios.au.dk>
Emne: Balladen om de 50% forklaret variation i Bennetsens analyse skyldes at R2 i excel med “fixed origin” beregnes forkert
Kære alle,
Siden vildtforsker Bennetsen (og dermed også vi andre) har været så plaget af at excel giver en R2 på 50%, hvor samme analyse i SAS giver 92%, har jeg lige foretaget en yderst grundig undersøgelse på Google over dette emne.
Søgeordene var ” r2 value linear regression excel fixed origin”. Efter lang tids søgen (2-4 minutter – jeg skulle lige tjekke et par andre hits først) fandt jeg følgende kilde, som forklarer det jeg selv antog i går, nemlig at der er en programmeringsfejl i excel 2010 og yngre version, som giver forkert værdi når R2 beregnes ved lineær regression fikseret skæring i (0,0):
http://www.quepublishing.com/articles/article.aspx?p=2019170
Postyret omkring de 50% burde dermed kunne afblæses.
Jeg håber I vil være behjælpelige med at videreformidle dette budskab til de af jeres medlemmer, der måtte have interesse deri.
Med venlig hilsen
Peter
Bilag 5:
http://www.quepublishing.com/articles/article.aspx?p=2019170 (Mine fremhævelser og tilføjelser):
Setting the const argument to FALSE can easily have major implications for the nature of the results that LINEST() returns. And there is a real question of whether the const argument is a useful option at all. In fact, the question is not limited to LINEST() and Excel. It extends to the whole area of regression analysis, regardless of the platform used to carry out the regression.
Some credible practitioners believe that it’s important to force the constant to zero in certain situations, usually in the context of regression discontinuity designs.
Others, including myself, believe that if setting the constant to zero appears to be a useful and informative option, then linear regression itself is often the wrong model for the data.
Note
I should mention that it’s easy to reach the conclusion that forcing the constant to zero results in a more accurate outcome. That belief is based on a higher value for R2, and an F ratio that argues more strongly for rejecting a null hypothesis of no relationship between the Y variable and the X composite. This (spurious)
conclusion is based on a misunderstanding of what happens mathematically when the constant is forced to zero.
The Constant and the Deviations
Of course, the problem is due to the fact that in omitting the constant, we are redefining what’s meant by the term “sum of squares.” As a result, we’re dismembering the meaning of the R2.
When you include the constant, the deviations are the differences between the observed values and their mean—that’s what “least squares” is all about. When you omit the constant, the deviations are the differences between the observed values and zero—that’s what “regression without the constant” is all about.
If the predicted values happen to be generally farther from zero than from their own mean, then the sum of squares regression will be inflated as compared to regression with the constant. In that case, the R2 will tend to be greater without the constant in the regression equation than it is with the constant.
Bilag 6:
From: Mads Flinterup [mailto:mf@jaegerne.dk]
Sent: 16 December 2016 12:15
To: Henrik Meilby
Cc: Peter Sunde; Niels Søndergaard
Subject: S1921 har behov for hjælp!
Importance: High
Hej Henrik
Du har før redet min røv ud i databehandling – det var dog i min tid som simpel forststuderende med studie nummeret S1921…
Sagen er den, at min karriere har ført mig derhen, at jeg bruger min tid på forvaltning af naturressourcen kronvildt, og i den forbindelse har et meget tæt samarbejde med Peter Sunde fra Institut for Bioscience på Aarhus Universitet. Kronvildtforvaltning er et meget hedt jagtpolitisk tema, og i den slags debatter finder modstanderne altid et eller andet vilkår at hænge deres hat på.
Jeg skal sparer dig for en helt masse udenoms ævl og komme til mit egentlig ærinde.
Over fem jagtsæsonen indsamlede vi kæber fra nedlagte krondyr her på Djursland. På baggrund heraf har man blandt andet beregnet hjortes overlevelsesrater m.v. for at kunne dette har vi naturligvis skulle fastsætte alderen på de nedlagte krondyr. I krondyrets første og andet leveår lader alderen sig entydigt definere via tandskiftet. Men når dyret har sit blivende tandsæt fra det tredje leve år er der behov for en anden metode. Her støtter man sig til, at der i tandroden dannes ringe i tandcementen. Når tandroden behandles rigtig og snittet læses kan man derved aflæse dyrets alder lige som årringe i et træ.
Nu er den øvelses på træer i nogle tilfælde nemmere end andre, og usikkerheden på et spredtporret træ markant større end på et hurtigt voksende ringporret træ.
For at belyse netop denne usikkerhed har man indsamlet kæber fra 38 krondyr med kendt alder og fået deres alder bestemt via tandsnitsmetoden. Rådata fremgår i vedhæftede excel fil som rådata til figur 5. Bemærk at tandsnitbedømmelse kun har været mulig på 37.
Striden står nu om hvor stor usikkerhed, der er på tandsnitsmetodens bestemmelse. Når lægmænd laver liniær regression på data for de en r2 væsentlig forskellig fra den videnskabelige r2. (Der er en fejl i Excel…)
Kan jeg bede dig om at køre en simpel lineær regression på de 37 observationer hvor du finder r2 på x=y og altså linjen skærer i 0,0?
I det øvelsen skulle have et omfang der kræve fakturering er jeg klar på det. Det slut produkt jeg skal stå tilbage med er at jeg kan sige at kompetent databehandler ved IFRO KU har gennem regnet data og fået r2 til…
Med venlig hilsen
Mads Flinterup
Schweiss og Hjortevildtkonsulent / Sektionsleder
Bilag 7:
From: Mads Flinterup [mailto:mf@jaegerne.dk]
Sent: 16. december 2016 19:10
To: Henrik Meilby
Cc: Peter Sunde; Niels Søndergaard
Subject: Sv: S1921 har behov for hjælp!
Fantastisk Henrik!
Indser pludseligt hvor meget jeg har savnet vores snakke.
Sagen er jo entydigt her at statistik er et sprog, der kun tales flydende af få personer. Personlig mestre jeg det efterhånden kun gebrokkent. (Imodsætning til andre sprog tror jeg faktisk ikke performance stiger med alkohol indtag!). Jeg er i øvrigt ikke i tvivl om at du og Peter ville kunne gennemføre en samtale på dansk, hvor jeg ville fatte uendeligt lidt.
Uanset hvad så står diskussionen om r2 værdien. De øvrige værdier er naturligvis også beregnet og anvendt i rapportens konklusioner.
Data er valideret. Dyrerne har kendt alder og metoden har givet svaret. Hvorfor NINA (de er bestem i Norge) har valgt at sige én er ubrugelig og så give to åbenlyse flyers ved jeg ikke. Men NINA har ikke kendt den rigtige alder, og de har altså opgive den alder de har bestemt ved snit.
Så analysen skal altså gøres på det komplette set.
Du er velkommen til at give mig alle de statistiske betragtninger du vil, men jeg skal kun bruger r2 ved skæring i 0,0 i forhold til Jægerforbundets pressemeddelelse i samarbejde med DCE i næste uge.
God weekend
Mads
