Afsnit 2.7: Øvelse 2: Geologi

Denne uges øvelse vedrører inferens i binomialmodellen og i poissonmodellen. I skal lave test i binomialmodellen og konfidensinterval i både binomialmodellen og poissonmodellen. I den første opgave skal I træne i at afgøre, om binomialmodellen eller poissonmodellen skal bruges til at beskrive data. Husk altid i opgaverne at opskrive den statistiske model til beskrivelse af data.
Opgaverne 2.1-2.5 skal være forberedt hjemmefra og gennemgås ved tavlen til øvelserne. Efter øvelsen skal der afleveres en rapport over opgave 2.6.

Opgave 2.1: Identificere fordeling

Betragt følgende fem datasæt. Forklar i hvert enkelt tilfælde, hvorfor du vil beskrive data med enten en binomialfordeling eller en poissonfordeling. Opskriv den statistiske model for hver observation.
  1. I et snit af en kalksten er der lagt et gitter ned over, og blandt 100 gitterpunkter er der i 28 fundet en bestemt fossil. Hvilken fordeling vil du bruge til at beskrive observationen ?
  2. I en jordprøve på 8 kilo er der fundet 17 guldkorn. Hvilken fordeling vil du bruge til at beskrive observationen ?
  3. På et geologihold på Aarhus Universitet er der startet 23 studerende, men efter et år har 4 valgt at skifte studie til kemi. Hvilken fordeling vil du bruge til at beskrive observationen ?
  4. I New Zealand har der i en periode på 30 år været 7 jordskælv med en styrke over 6 på Richterskalaen. Hvilken fordeling vil du bruge til at beskrive observationen ?

Opgave 2.2: Test i binomialmodellen

I Opgave 1.5 er omtalt en opdeling af jordskælv i forhold til tidevandsfasen. Af de 266 jordskælv er der 117, der falder under den opadgående fase af tidevandet og 149 i den nedadgående fase.
  1. Start med den tænkte situation, at der ud af 266 jordskælv er observeret 133 i den opadgående tidevandsfase. Opstil binomialmodellen til beskrivelse af denne tænkte situation. Forklar, at -værdien er lig med 1 for test af hypotesen, at sandsynlighedsparameteren i binomialmodellen har værdien
  2. Betragt nu den faktiske datasituation med 117 observerede jordskælv i den opdgående tidevandsfase. Opstil binomialmodellen for dette eksperiment og find -værdien for test af hypotesen, at sandsynlighedsparameteren i binomialmodellen har værdien (benyt eventuelt R-koden nederst i afsnit 1.3). Hvordan passer resultatet med dine simulationer i Opgave 1.5?
I Afsnit 1.1 omtalte jeg et af Mendels eksperimenter. Mendel undersøgte mange egenskaber ved ærteplanten. I tabellen nedenfor er et udsnit af Mendels eksperimenter. For alle de viste data er Mendels hypotese, at sandsynligheden for en recessiv er
  1. Betragt binomialmodellen for hver række i tabellen, og find de manglende -værdier for test af hypotesen, at sandsynligheden for recessiv er Hvor mange af de syv -værdier er under 0.5?
Hvis alle hypoteser er sande, forventer vi, at -værdien cirka halvdelen af gangene vil være under 0.5 og halvdelen af gangene være over 0.5. Mendels eksperimenter har været kritiseret for, at -værdierne generelt er for høje. I kan finde en diskussion af dette i artiklen Are Mendel's Data Reliable?.

Opgave 2.3: Konfidensinterval i binomialmodellen

Betragt data fra spørgsmål (b) i Opgave 2.2.
  1. Opstil en statistisk model for data, og lav et 95%-konfidensinterval for sandsynligheden for, at et jordskælv falder i den opadgående tidevandsfase (benyt eventuelt R-koden i Eksempel 2.2.3).

Opgave 2.4: Konfidensinterval i poissonmodellen

Antallet af jordskælv i et bestemt område og med en styrke i et givet interval beskrives ofte med en poissonfordeling. Et eksempel er artiklen A Poisson model for earthquake frequency uncertainties in seismic hazard analysis. I artiklen betragtes blandt andet jordskælv i New Zealand. Information om disse kan findes på nettet under adressen info.geonet.org.nz. Data i tabellen nedenfor viser antallet for tre styrkeintervaller og for perioden 1930-2015. Styrken er på Richterskalaen, som er en logaritmisk skala. Hvis styrken stiger med 1, stiger den samlede energi i jordskælvet med Gutenberg-Richter loven for jordskælv angiver forholdet mellem antallet af jordskælv af forskellig styrke. I tabellen er dette forhold angivet for de tre styrkeintervaller (med "b-value" i Gutenberg-Richter loven sat til 1).
  1. Betragt styrkeintervallet 6.0-6.3. Opskriv poissonmodellen for data, idet raten i poissonmodellen er det forventede antal jorkskælv per år for det betragtede styrkeinterval. Lav et 95%-konfidensinterval for raten (benyt eventuelt R-koden i Eksempel 2.5.4).
  2. Betragt styrkeintervallet 6.3-6.6. Hvis Gutenberg-Richter loven holder, er det ifølge tabellen ovenfor relevant at skrive raten per år som hvor så tolkes som raten per år for styrkeintervallet 6.0-6.3. Lav et 95%-konfidensinterval for baseret på data i styrkeintervallet 6.3-6.6.
    Hvad er din umiddelbare vurdering: er der samme værdi af raten her som i det foregående spørgsmål?
  3. Hvis vi tror på Gutenberg-Richter loven, kan vi kombinere data for alle tre styrkeintervaller i tabellen ovenfor for at lave et fælles konfidensinterval for raten per år i styrkeintervallet 6.0-6.3. Det samlede antal jordskælv er og raten per år er Opskriv poissonmodellen for det samlede antal jordskælv, og lav et 95%-konfidensinterval for

Opgave 2.5: Fortolkning af konfidensinterval

Denne opgave går ud på at eftervise fortolkningen af et konfidensinterval for poissonfordelte data ved hjælp af simulationer. For at løse opgaven skal I bruge koden i kodevinduet nedenfor. I R kan man simulere observationer fra en poissonfordeling med kommandoen rpois, Således vil kommandoen rpois(10,3) give 10 observationer fra en poissonfordeling med rate (middelværdi) lig med 3.
Til hjælp i opgaven er der i kodevinduet nedenfor angivet en kode, der simulerer observationer fra en -fordeling, og for hver simuleret værdi beregner det tilhørende 95%-konfidensinterval. I beregningen bruges funktionen mean, som beregner gennemsnit af data i en vektor.
Først følger en række spørgsmål til forståelse af koden.
  1. Hvilken poissonfordeling simuleres der fra i kodevinduet?
  2. Hvor mange observationer bliver der simuleret?
  3. Hvad beregnes i lower?
  4. Hvad repræsenterer den sidste værdi der udskrives?
  5. Prøv at køre koden et par gange og se variationen i resultatet.
Nu følger en række spørgsmål til belysning af et konfidensinterval.
  1. Simuler for lig med henholdsvis og en enkelt værdi (nSim sættes til 1) fra en poissonfordeling med middelværdi med raten
    Lav en tabel med fem søjler med henholdsvis værdien af den simulerede værdi fra poissonfordelingen, nedre og øvre grænse i konfidensintervallet og længden af konfidensintervallet. Tabellen har fire rækker svarende til de fire værdier af
    Kan du formulere en grov regel for, hvad der sker med længden af konfidensintervallet, når tælletiden firedobles?
  2. Betragt igen de fire tilfælde med og Benyt for hver værdi af koden ovenfor med nSim=100000.
    Lav en tabel med tre søjler med henholdsvis værdien af andel af de simulerede konfidensintervaller der indeholder den sande værdi af og middellængden af de nSim konfidensintervaller. Tabellen har fire rækker svarende til de fire værdier af
    Kan du genfinde din grove regel fra foregående spørgsmål, når du kigger på middellængden af konfidensintervallerne?
  3. I din tabel vil du se, at andelen af intervaller indeholdende den sande værdi af parameteren næsten ikke afhænger af Kan du forklare, hvorfor dette er tilfældet?

Opgave 2.6: Afleveringsopgave

Ved restaurering og vedligehold af vandløb i bymæssige omgivelser vil man gerne have viden om, hvor stabilt det opbyggede system er. For eksempel kan man lave en bund af store sten et sted, hvor vandløbet har et fald, og vil så gerne vide, om stenene bliver liggende eller flytter sig i perioder med stor vandføring. Dette kan man gøre ved at mærke nogle sten. Et sådant studie er beskrevet i Assessing the performance of in-stream restoration projects using radio frequncy identification (RFID) transponders. I artiklen foregår mærkningen ved, at der bores hul i stenen, og en sender indsættes i hullet. På denne måde kan man individuelt mærke sten. I artiklen omtales et eksperiment vedrørende restaurering af Wilket Creek, som er et sideløb til Don River i Toronto. I efteråret 2013 er der lagt mærkede sten ud to steder, henholdsvis det restaurerede område og et kontrolområde uden restaurering. I foråret 2014 har man fundet, hvor mange af de mærkede sten der har flyttet sig. Eksperimentet er gentaget henover august måned 2014. Resultaterne fremgår af tabellen nedenfor. I den første periode har vandstanden en enkelt gang været oppe på det maksimale niveau, som restaureringen er designet til.
Vi er interesseret i at lave skøn og konfidensinterval for den andel (populationsfrekvensen) af alle stenene i vandløbet, der potentielt vil flytte sig ved store vandføringer.
  1. Betragt data fra periode 1 og restaureringsområdet (prøvenummer 1). Angiv en statistisk model til beskrivelse af data, og beregn et 95%-konfidensinterval for andelen af sten, der flyttes.
  2. Lav en tabel, hvor skøn og konfidensinterval for andelen af sten, der flyttes, angives for alle fire prøvenumre i ovenstående tabel.
  3. Til sidst skal du lave en grafisk illustration af dine resultater. Lav en figur, hvor skøn over andel af sten, der flyttes, afsættes mod prøvenummer med kommandoen
    plot(proevenr,skoen,ylim=c(0,1))
    hvor proevenr er vektoren c(1,2,3,4), og skoen er en vektor med de fire skøn over andel af sten, der flyttes. Undersøg, hvad tilføjelsen ylim=c(0,1) betyder ved at prøve at ændre tallene 0 og 1.
    Du skal slutte af med at indtegne de fundne konfidensintervaller i figuren i form af lodrette linjestykker. Der er ikke nogen oplagt funktion i R til at gøre dette, og jeg har derfor selv kodet en funktion errorbar. Se punktet Egne funktioner i R i afsnit 1.2 med hensyn til adgang til denne funktion. Når funktionen er hentet ind i R, kan du indtegne konfidensintervallerne i figuren med kommandoen
    errorbar(proevenr,skoen,nedre,oevre)
    hvor nedre er en vektor med nedre endepunkter for de fire konfidensintervaller og oevre en vektor med de øvre endepunkter.
Bemærkning: Konfidensintervallet i en binomialmodel kan udregens i R med kommandoen prop.test(x,n,correct="FALSE")conf.int. Prøv denne funktion med en af dine udregninger ovenfor.

ForegåendeNæste