Formålet med data i denne opgave er at prædiktere, hvor god en kemisk forbindelse er til at bremse BTK (Brutons Tyrosine Kinase). Respons er $\text{tIC}_{50}=-\log_{10}(\text{IC}_{50})$, hvor $\text{IC}_{50}$ er half maximal inhibitory concentration, det vil sige den mængde, der er nødvendig for at nedsætte BTK aktiviteten til det halve ($\text{IC}_{50}$ skal angives i mol/L før logaritmen beregnes). Bemærk, at en stor værdi af $\text{tIC}_{50}$ betyder større evne til at bremse BTK. Hvis man ud fra simple beskrivelser af den kemiske forbindelse kan prædiktere $\text{tIC}_{50}$, kan man benytte dette til at designe nye forbindelser eller til at løbe gennem et bibliotek af kemiske forbindelse for at finde gode kandidater til at bremse BTK. Data i denne opgave er en delmængde af data fra Computational regression and prediction analysis on a dataset of 52 Pyrrolo[2,3-b]Pyridines and Pyrimidines Rheumatoid Arthritis inhibitors. Vi vil betragte fire variable til beskrivelse af $\text{tIC}_{50}$: molekylevægten af den kemiske forbindelse (MW), logaritmen til partition coefficient (logP), antal hydrogenbindinger i acceptorgruppe (HBA) og måling foretaget ved røngtenspektroskopi (KA1). Forfatterne identificerer to datapunkter som outliers, og data, I skal betragte, indeholder ikke disse. Data findes i filen tICdata.csv, som har 5 søjler i rækkefølgen MW, logP, HBA, KA1 og tIC50. Der er kun 48 datarækker i filen, idet jeg har taget to sæt ud, som I skal bruge til sidst i opgaven til prædiktion.

1. Indlæs data som en datatabel.
   Opskriv den fulde regressionsmodel, hvor middelværdien af tIC50 afhænger lineært af de fire forklarende variable.
2. Lav et qqplot af residualerne i denne model. Angiv Resultat nummer for hvordan et qqplot vurderes
3. Reducer den fulde multiple regressionsmodel ved successivt at fjerne led i modellen (backward selektion). Lav en tabel, som for hvert trin angiver den variabel der fjernes i modellen og den tilhørende $p$-værdi.
4. Opskriv slutmodellen ved backward selektionsproceduren, og lav et $F$-test for reduktion fra den fulde model til slutmodellen.
5. Lav en figur, hvor residualerne fra slutmodellen afsættes mod KA1, og med nullinjen indsat. (Overvej også, hvilke andre figurer du burde lave.) Idet du starter fra venstre side og går mod højre, hvor mange punkter ligger under nullinjen inden det første punkt, der ligger over nullinjen ?
6. Lav en figur, hvor tIC50 afsættes mod de forventede værdier, og indsæt identitetslinjen i denne figur.
   Synes du, at slutmodellen vil være god til at prædiktere tIC50 ?
7. Angiv 95%-konfidensintervaller for de to regressionsparametre i din slutmodel. Hvilke kemiske forbindelser skal man lede efter, hvis man gerne vil have en forbindelse, der er god til at bremse BTK ?
8. Datasætet indeholder to observationer, der ikke er medtaget i data, som I benyttede ovenfor. Data for de to observationer er
   $$\begin{array}{ccccc} \hline \text{MW}& \text{logP}& \text{HBA}& \text{KA1} & \text{tIC50} \\ \hline 573.67 & 5.4525 & 5 & 31.5918 & -1.5224 \\ 578.7 & 3.4987 & 5 & 31.5918 & -0.8451 \\ \hline \end{array}$$ Lav et 95%-prædiktionsinterval for tIC50 for hver af de to prøver ovenfor. Ligger den observerede værdi i prædiktionsintervallet ?
   Som nævnt i indledningen til denne opgave, blev to prøver fjernet fra datasættet, idet disse blev betragtet som outliers. Data for disse to er
   $$\begin{array}{ccccc} \hline \text{MW}& \text{logP}& \text{HBA}&\text{KA1} & \text{tIC50} \\ \hline 497.57 & 3.2129 & 5 & 27.1837 & -3.1326 \\ 580.72 & 3.7191 & 5 & 31.5918 & -2.7589 \\ \hline \end{array}$$ Gentag ovenstående beregning af prædiktionsinterval for disse to prøver.

Luk igen

I denne øvelse skal I se på NIR-spektre af småkagedej for at vurdere, om disse kan bruges til at prædiktere indholdet af vand i dejen. Data er hentet fra ppls pakken i R (programpakke til statistikberegninger), og er analyseret i de to artikler Application of near-infrared reflectance spectroscopy to compositional analysis of biscuits and biscuit dough, og Bayesian Wavelet Regression on Curves with Applications to a Spectroscopic Calibration Problem. Der er 40 prøver med NIR-spektre målt ved 700 bølgelængder (fra 1100 til 2498 nanometer). Data ligger i filen DejTraen.txt, som har 701 søjler (uden søjleoverskrifter), hvor den sidste søjle indeholder mængden af vand. Hver række i filen svarer til en prøve.

1. Indlæs data som en $40\times 701$ matriks Traen med brug af indlæsningskommandoen np.loadtxt.
   Placer NIR-spektrene i en $40\times 700$ matrix $T$ og mængden af vand i vektoren vand (python: T=Traen[:,0:700]).
   Angiv, hvor mange dejprøver der er med et vandindhold under 13.
2. Tabellen nedenfor viser resultatet af et kald til forward samt et kald til cvForward.
   $$\begin{array}{lcccccc} \hline \text{Antal variable} & 1 & 2 & 3 & 4 & 5 & 6 \\ \hline \text{Variabelnummer} & 425 & 588 & 641 & 242 & 11 & 10 \\ \text{Spredningsskøn} & 0.924 & 0.543 & 0.525 & 0.486 & 0.441 & 0.424 \\ \text{P-værdi} & 0.0000 & 0.0000 & 0.0656 & 0.0124 & 0.0059 & 0.0617 \\ \hline \text{CV-prædiktionsspredning} & 1.248 & 0.870 & 0.975 & 0.832 & 0.870 & 0.879 \\ \hline \end{array}$$ Illustrer forward selektionsprocessen med en figur som i afsnit 9.4 baseret på data i tabellen ovenfor.
   Hvor mange variable, fundet ved forward selektion, vil du medtage i en model til beskrivelse af vandindholdet ?
3. Du skal nu bruge ridge regression til at beskrive vandmængden ud fra NIR-spektrene. Lav en tabel med krydsvalideringsspredning $s_{\text{cv}}(\lambda)$ for $\lambda$-værdierne $0.001,0.01,0.1,1,10,100,200,400,1000$. Lav også en figur, hvor $s_{\text{cv}}(\lambda)$ afsættes mod $\log(\lambda)$.
   Vælg en værdi $\hat\lambda$ af $\lambda$, som du vil bruge til at beskrive data med. Lav en figur med den observerede vandmængde afsat mod den forventede vandmængde, og indsæt identitetslinjen i denne figur. Du kan eventuelt lade dig inspirere af koden i den skjulte kode Ridge regression i python.
4. Hvilken model foretrækker du, modellen fra forward selektion eller modellen fra ridge regression ?

Luk igen

Denne opgave er en fortsættelse af den foregående opgave. Udover det oprindelige datsæt med 40 dejprøver er der også et nyt datasæt med 32 prøver. Vi kan bruge dette datasæt til et uafhængigt tjek af, hvor god modellen er til at prædiktere. Det nye datasæt er i filen DejTest.txt, som har samme struktur som DejTraen.txt (men altså 32 rækker i stedet for 40 rækker).

1. Indlæs som i foregående opgave data fra filen DejTraen.txt og dan en matriks $T$ med spektrene og en vektor vand med indholdet af vand for de 40 prøver.
   Indlæs dernæst data fra filen DejTest.txt og dan en matrik Ttest med spektrene og en vektor vandtest med indholdet af vand for de 32 nye prøver.
2. Betragt først forward selektionsmodellen fra foregående opgave, hvor der medtages 4 variable. Estimer modellen baseret på $T$ og vand, og lav prædikterede værdier for prøverne med spektre i Ttest. Dette kan gennemføres under brug af koden i det skjulte punkt Bruge slutmodellen til prædiktion.
   Udregn testspredningen som kvadratroden af den gennemsnitlige værdi af den kvadrerede afstand mellem prædiktion og den sande værdi, det vil sige
   $$s_{\text{test}}=\sqrt{\frac{1}{32} \sum_i (\text{vandtest}_i-\text{prædiktion}_i)^2}.$$ Overrasker værdien dig ?
3. Betragt nu ridge regression, hvor regulariseringsparameteren $\lambda$ vælges til $\lambda=3.$ Benyt funktionen ridge til at finde $\hat\alpha(\lambda),$ $\hat\beta(\lambda)$ og funktionen cvRidge til at finde prædiktionsspredningen $s_{\text{cv}}(3)$ baseret på spektrene i $T$ og respons i vand.
   Angiv den første koordinat i vektoren $\hat\beta(\lambda).$
4. Udregn de prædikterede værdier for de nye data med spektre i Ttest (se koden i det skjulte punkt Ridge prædiktion i python)
   Lav en figur hvor de observerede værdier for de nye data afsættes mod de prædikterede værdier. Indsæt identitetslinjen i figuren.
   Hvad er den største abolutte prædiktionsfejl, og hvor optræder den ?
5. Udregn testspredningen som kvadratroden af den gennemsnitlige værdi af den kvadrerede afstand mellem prædiktion og den sande værdi.
   Sammenlign denne værdi med cross-validation prædiktionsspredningen $s_{\text{cv}}(3).$ Synes du ridge regressionsmodellen er god til at prædiktere nye dejprøver ?

Luk igen

I opgave 5.8 så I på en regressionsmodel for absorbansratioens afhængighed af proteinmængden BSA, og brugte sammenhængen som en kalibreringskurve. I opgaven var der givet en enkelt observation af absorbansratioen for hver af seks værdier af BSA. I artiklen Linearization of the Bradford Protein Assay, som data stammer fra, er hver absorbansratio faktisk gennemsnit af tre uafhængige målinger. Alle data fra artiklen er gengivet i tabel nedenfor og findes i filen TreGentagBradford.csv. Filen har tre søjler, hvor første søjle indeholder BSA, den anden søjle er en variabel gruppe, der koder for BSA-gruppe (med værdierne G1,G2,G3,G4,G5,G6), og den tredje søjle indeholder absorbansratio. $$\begin{array}{c|ccc|ccc|ccc}\hline && \text{590nm} &&& \text{450nm} &&& \text{F} & \\ \text{BSA} & \text{BrønA} & \text{BrønB} & \text{BrønC} & \text{BrønA} & \text{BrønB} & \text{BrønC} & \text{BrønA} & \text{BrønB} & \text{BrønC} \\ \hline 0 & 0.290 & 0.298 & 0.306 & 0.629 & 0.654 & 0.672 & 0.461 & 0.456 & 0.455 \\ 1 & 0.491 & 0.457 & 0.462 & 0.600 & 0.593 & 0.590 & 0.818 & 0.771 & 0.783 \\ 2 & 0.580 & 0.591 & 0.601 & 0.516 & 0.534 & 0.539 & 1.124 & 1.107 & 1.115 \\ 3 & 0.683 & 0.717 & 0.751 & 0.483 & 0.495 & 0.495 & 1.414 & 1.448 & 1.517 \\ 4 & 0.802 & 0.853 & 0.840 & 0.446 & 0.459 & 0.444 & 1.798 & 1.858 & 1.892 \\ 5 & 1.226 & 0.965 & 0.971 & 0.677 & 0.394 & 0.403 & 1.811 & 2.449 & 2.409 \\ \hline \end{array}$$ Et relevant spørgsmål for data af denne type er, om man skal tage gennemsnit for hver bølgelændge af de tre målinger, før man beregner forhold, eller om man skal lave forholdet for hver brønd og tage gennemsnit af disse. Formodentligt giver det en meget lille forskel, men da respons er lineær for forholdet, vil det være bedst at betragte forholdet for hver brønd og så tage gennemsnit. Et andet spørgsmål I skal overveje nedenfor er, om man skal basere kalibreringskurven på gennemsnit af de tre forhold for hver proteinmængde, eller om man skal basere beregningerne på alle 18 datapunkter.

1. Opstil normalfordelingsmodellen, hvor hver BSA-gruppe har sin egen middelværdi af absorbansratio, og alle grupperne har den samme varians.
   Analyser modellen i python, og angiv skøn over spredningen i modellen.
2. Opstil regressionsmodellen, hvor middelværdien af absorbansratio afhænger lineært af BSA.
   Analyser modellen i python, og angiv skøn over spredningen i modellen.
   Hvis man kun betragter data bestående af 6 proteinmængder og 6 absorbansratioer, som hver er gennemsnit af tre målinger, vil en regressionsmodel give et skøn over spredningen på 0.02665. Kan du forstå dette i sammenligning med det skøn, du fik ovenfor (svaret kan godt være "nej") ?
3. Forklar, at regressionsmodellen er en undermodel af modellen, hvor hver BSA-gruppe har sin egen middelværdi af absorbansratio (den ensidede variansanalysemodel).
   Lav et test for reduktion fra den ensidede variansanalysemodel til regressionsmodellen. Angiv den $F$-fordeling, der bruges. Er $p$-værdien stor nok, til at du vil acceptere regressionsmodellen ?
4. Beregn et 95%-konfidensinterval for hældningen i den lineære sammenhæng mellem middelværdien af Absorbansratio og BSA.
   Hvis man igen betragter regressionsmodellen baseret på de 6 gennemsnit, bliver konfidensintervallet $[0.335,0.371].$ Sammenlign dette med det konfidensinterval du lige har fundet.
5. Betragt en ny prøve med ukendt indhold af proteinmængden BSA, hvor absorbansratio er målt tre gange med værdierne 0.979, 0.977, 0.984. Angiv et 95%-konfidensinterval for det ukendte indhold af BSA.
   Hvis man igen betragter regressionsmodellen baseret på de 6 gennemsnit, og en ny måling som er 0.980 (gennemsnittet af de tre målinger ovenfor), bliver konfidensintervallet $[1.31,1.78].$
   Ved fremtidig brug af kalibreringskurven til bestemmelse af proteinindholdet vil du da vælge at bruge kurven og spredningskøn baseret på alle 18 datapunkter, eller vil du bruge kurven og spredningsskøn baseret på de 6 gennemsnit ?

Luk igen

I denne sidste opgave er det meningen, at I hovedsageligt skal bruge metoderne fra de 6 første øvelser (undtagelsen er det sidste spørgsmål). I artiklen Use of Multiple Linear Regression Models for Setting Water Quality Criteria for Copper: A Complementary Approach to the Biotic Ligand Model betragtes forskellige modeller til at beskrive giftigheden af kobber i vandløb. Giftigheden, under givne forhold, findes ved at lave en række forsøg med forskellige koncentrationer af kobber i det vand hvor vandlopper (Daphnia magna) opholder sig i et givet tidsrum. Som beskrivelse af giftigheden bruges den koncentration af kobber under hvilken 50 procent af vandlopperne dør. Denne kaldes acute EC50 (Half maximal effective concentration). Som modeller for EC50-koncentrationen har man på den ene side en simplificeret model, hvor man kun bruger vandets hårdhed som forklarende variabel i en regressionsmodel. På den anden side har man en kompliceret multipel regressionsmodel med 10 forklarende variable, som ifølge forfatterne ikke har vundet indpas i vandkvalitetsprogrammer. Forfatterne betragter i stedet en mere simpel multipel regressionsmodel med tre forklarende variable, henholdsvis DOC (dissolved organic carbon, mg/L), hårdhed (mg/L) og pH. I opgaven her skal I betragte data for EC50-koncentrationen af kobber ($\mu$g/L Cu). Den bedste skala at betragte data på er ved at logaritmetransformere værdierne (dog ikke pH, som allerede er på en log-skala). Filen Dmagna.csv indeholder variablene DOC, hardness, pH og EC50. Der er data for 302 målinger. Lad logEC50 være logaritmen til acute EC50, og for nemhed i notationen lad $D$ være logaritmen til DOC og $H$ være logaritmen til hardness. En analyse af den multiple regression af logEC50 på $D$, $H$ og pH viser, at modellen kan reduceres til $$\mathrm{LogEC50}\sim N(\alpha+D+\frac{1}{2}H+\mathit{pH},\sigma^2).$$ Du kan eventuelt køre en multipel regressionsanalyse og se, at parametertabellen ikke strider mod $\beta_D=1,$ $\beta_H=\frac{1}{2}$ og $\beta_{\text{pH}}=1.$ Denne analyse giver anledning til at indføre et indeks, her kaldet giftIndex, på formen $$\text{giftIndex}=\text{logEC50}-D- \frac{1}{2}H-\mathit{pH}.$$ Fra data ses, at dette indeks svinger omkring værdien $-6.86.$ Opgaven her går ud på at lave nogle (lidt usædvanlige) undersøgelser for at se, om den multiple regressionsmodel kan forbedres ved at inkludere flere variable. Specifikt skal I se, om produktet $H\cdot\mathit{pH}$ indeholder information om giftindex ved at dele data op i grupper efter denne variabel.

1. Indlæs data og dan variablene logEC50, $D$, $H$ og pH. Dan dernæst giftIndex.
   Tjek at gennemsnit af målingerne af giftIndeks er $-6.86.$
2. Konstruer en kategorisk variabel gruppe, der deler op i tre grupper, $A$, $B$ og $C$, efter om $H\cdot\mathit{pH}$ er mindre end 25, er mellem 25 og 35 eller over 35. I python (pandas) kan dette gøres med kommandoen
   $$\text{gruppe=pd.cut(H*pH,bins=[0,25,35,100],labels=['A','B','C'])}$$
   Dan nu to datasæt, giftA og giftC, med værdierne af giftIndex hørende til gruppe $A$ og til gruppe $C$. Der er 60 observationer i gruppe $A$ og 96 i gruppe $C$.
3. For de 60 prøver i gruppe $A$ er der 14, der har en værdi af giftIndex under $-7$ (konstater, at dette er rigtigt). Opstil en statistisk model til beskrivelse af observationen 14, og lav et 95%-konfidensinterval for sandsynligheden for, at værdien af giftIndex er under $-7$ i gruppe $A$.
4. For de 96 prøver i gruppe $C$ er der 29, der har en værdi af giftIndex under $-7.$ Undersøg, om der er samme frekvens af prøver med en værdi af giftIndex under $-7$ i de to grupper $A$ og $C$.
5. Opstil en statistisk model til beskrivelse af værdierne af giftIndex i giftA og giftC.
   Lav et test for hypotesen, at der er samme middelværdi af giftIndex i de to grupper $A$ og $C$.
6. Betragt nu giftIndex for alle 3 grupper dannet ud fra faktoren gruppe. Opstil en statistisk model for data, og undersøg først, om der er samme varians for de tre grupper, og dernæst, om der er samme middelværdi for de tre grupper.
7. Ovenstående analyse viser, at der er information i produktet $H\cdot\mathit{pH}$, der kan bruges til at beskrive giftIndex, og dermed kan bruges til at beskrive logEC50. I dette delspørgsmål skal I analysere en multipel regressionsmodel til beskrivelse af logEC50, hvor I som forklarende variable ud over $D$, $H$ og pH også bruger alle produkter af disse, det vil sige, variablene
   $$\begin{aligned} & \mathit{DD}=D*D,\enspace \mathit{HH}=H*H,\enspace \mathit{ppH}=\mathit{pH}*\mathit{pH}, \\ & \mathit{DH}=D*H,\enspace \mathit{DpH}=D*\mathit{pH},\enspace \mathit{HpH}=H*\mathit{pH}. \end{aligned}$$ Opskriv den multiple regressionsmodel, hvor alle 9 forklarende variable inddrages.
   Reducer modellen ved brug af backward selektion.
8. Lav et test for reduktion fra startmodel til slutmodellen fremkommet ved backward selektion.
   Angiv et 95%-konfidensinterval for regressionskoefficienten hørende til variablen $H\cdot\mathit{pH}$ i slutmodellen.
   I din slutmodel, hvor meget stiger middelværdien af logEC50 hvis dissolved organic carbon stiger fra 2 til 10 ?
9. Overvej, hvilke figurer du ville lave til kontrol af slutmodellen (du skal ikke lave figurerne).
10. Sammenlign endelig skøn over spredning i slutmodellen med skøn over spredning i den multiple regressionsmodel hvor kun de tre forklarende variable $D$, $H$ og pH inddrages. Er det rimeligt at bruge denne sidste simple model ?

Luk igen

Betragt data fra opgave 7.1 og den multiple regressionsmodel for tIC50 med de to forklarende variable logP og KA1. Analyse af denne model giver $$\text{Cov}_s(\hat\alpha,\hat\beta_{\text{logP}})=-0.01813,\enspace \text{Cov}_s(\hat\alpha,\hat\beta_{\text{KA1}})=-0.01029,\enspace \text{Cov}_s(\hat\beta_{\text{logP}},\hat\beta_{\text{KA1}})=0.00007653.$$

1. Lav et approksimativt 95%-konfidensinterval for parameteren $\theta=\alpha+30\big(\beta_{\text{logP}}^2+\beta_{\text{KA1}}^2\big).$

Luk igen