1. No category

Vækstmodeller
Henrik L. Pedersen
Institut for Grundvidenskab og Miljø
Det Biovidenskabelige Fakultet
Københavns Universitet
16. september 2010
1
Matematiske modeller generelt
Matematiske modeller i naturvidenskab omhandler beskrivelse af virkelighedens naturvidenskabelige fænomener vha. matematik. Det kræver kendskab til et naturvidenskabeligt fagområde, til matematik og statistik og til it. Formålet med denne case er at give eksempler på
modellering inden for det biovidenskabelige område. Hovedeksemplet er dosis-respons kurver og særligt et fænomen, der går under betegnelsen hormesis. Som det vil fremgå kræver en
forståelse af dette kendskab til økotoksikologi, matematik (funktionsundersøgelse), statistik
(kurvetilpasning foruden metoder til parameterestimater, standardafvigelser og signifikans,
der tages op senere) og it (beregninger og visualisering). Vi vil bruge Maple som it-værktøj i
denne case, skønt statistikprogrammet R nok er mere naturligt at bruge (R introduceres senere
i uddannelsen). Fokus i denne case vil være at belyse hvordan en given model beskrives vha
matematiske værktøjer og visualiseres ved it.
En matematisk model skal på den ene side være så simpel, at den kan håndteres og på den
anden side være så kompliceret, at den giver et nogenlunde retvisende billede af virkeligheden.
En matematisk model kan (i varierende grad) være teoretisk funderet, dvs. baseret på veletablerede teorier inden for fx fysik eller kemi. Modeller kan også i højere grad være empirisk funderede – og her er det ikke altid muligt at begrunde modellens formler teoretisk, men i højere
grad at vælge på basis af observerede data. Dette er ofte tilfældet med matematiske modeller
inden for biologividenskab, fordi målingerne ofte vil være et resultat af en kombination af
mange (og komplicerede) processer.
Forbindelsen mellem på den ene side matematiske funktioner og på den anden side observerede data er først og fremmest kurvetilpasning. Vi beskriver derfor kort dette værktøj.
1
2
Kurvetilpasning
Igennem to forskellige punkter i XY-planen går netop en ret linje! Dette er utvivlsomt rigtigt,
men måske ikke så let at arbejde med, hvis man har en mængde punkter i planen som ser ud til
næsten, og kun næsten, at ligge på en ret linje. Hvordan afgør man hvilken ret linje der passer
bedst med punkterne? Svaret er: man benytter lineær regression!
For at beskrive, hvad denne metode går ud på skal vi have lidt notation. Vi tænker os, at der er
foretaget et antal (n) målinger (t1 , y1 ), . . . , (tn , yn ) og vi vil ud fra disse målinger bestemme en
ligning y = pt + q for den rette linje, der afviger mindst fra målingerne.
Hvordan skal vi beskrive afvigelsen mellem en linje og en række målinger? Dertil beregnes for
hver måling (ti , yi ) afvigelsen pti + q − yi mellem den målte værdi yi og den beregnede værdi
pti + q. Som et mål for denne afvigelse benyttes kvadratet ( pti + q − yi )2 og vi er interesserede
i at bestemme værdier for p og q sådan at summen af alle disse kvadratafvigelser
n
Q( p, q) =
∑ ( pti + q − yi )2
i =1
er så lille som muligt. Metoden illustreres i nedenstående eksempel. Se [4].
Eksempel 2.1 I et forsøg undersøges en mulig sammenhæng mellem luftfugtigheden L (målt i %) og
udklækningstiden U (målt i timer) for en type flueæg. Data findes i tabellen neden for.
L
U
100
16.6
94
17.4
88
18.3
82
18.2
76
19.6
70
20.2
64
21.7
58
22.0
52
22.7
46
23.2
Maple har funktionalitet der dækker lineær regression.
with(Statistics);
ler := [100, 94, 88, 82, 76, 70, 64, 58, 52, 46];
uer := [16.6, 17.4, 18.3, 18.2, 19.6, 20.2, 21.7, 22.0, 22.7, 23.2];
line := LinearFit([1, t], ler, uer, t);
luer := [seq([ler[k], uer[k]], k = 1 .. 10)];
plot(luer, style = point, symbol = asterisk, color = black,
labels = [Udklækningstid(L), Luftfugtighed(U)],
labeldirections = [horizontal, vertical]);
plu := %
plot(line, t = 40 .. 100); pline := %;
with(plots);
display([plu, pline]);
Køres disse linjer i Maple findes regressionslinjen til U = 29.29 − 0.13L, og det svarer til p = −0.13
og q = 29.29. Figur 1 viser resultatet af display([plu, pline]);
2
Figur 1: Målepunkter og linjen U = 29.29 − 0.13L
Som det ses beskrives data meget godt af en ret linje. Den centrale kommandolinje i Maple ovenfor er
kaldet LinearFit, der bevirker, at parametrene p og q beregnes. (Bemærk, at kaldet, der starter med
plot(luer af typesætningshensyn er delt to steder.)
Det viser sig, at de bedste værdier af p og q faktisk kan bestemmes som løsningen til to ligninger med to
ubekendte, som kan formuleres vha. matricer og vektorer. Vi vil dog ikke gå nærmere ind i disse detaljer. I
stedet for vil vi bruge Maple til at visualisere sammennhængen mellem forskellige værdier af parametrene
p og q og den tilsvarende afvigelse Q( p, q). Det kan gøres ved at plotte Q( p, q) som en funkton af to
variable:
with(plots):
Q := (p,q)->sum((p*ler[j]+q-uer[j])^2, j = 1 .. 10);
plot3d(Q(p,q), p=-0.5..0.5, q=10..50, axes=framed);
Resultatet af disse linjer ses i Figur 2.
I ovenstående eksempel lå data pænt omkring en ret linje. Ofte er man i den situation, hvor
man slet ikke forventer, at der kunne tænkes en lineær sammenhæng i data. Det kunne være
modeller for radioaktivt henfald, hvor ændringen af et stofs koncentration naturligt er proportional med stoffets koncentration. I Eksempel 2.2 (se [4]) giver vi et andet eksempel på en ikke
lineær samenhæng mellem pattedyrs kropsvægt og stofskifte.
3
Figur 2: Illustration af Q( p, q)
Eksempel 2.2 I tabellen neden for er angivet samhørende målinger af forskellige pattedyrs kropsvægt
(målt i kg) og stofskifte (målt i antal liter ilt forbrændt pr. time).
Dyr
Kropsvægt, K
Stofskifte, S
Vampyrflagermus
0.029
0.027
Næsebjørn
3.9
1.0
Ørkenræv
1.1
0.4
Hyæne
Kænguru
Jordsvin
7.0
2.2
33
5.8
48
6.0
Plottes disse målinger i et KS-koordinatsystem ses det, at data kun dårligt ville kunne beskrives ved en
lineær sammenhæng. Det viser sig, at en potenssammenhæng langt bedre beskriver data.
Vi vil nu bestemme parametrene a og b således at funktionen f (K ) = bK a kan siges at beskrive data
“bedst muligt”.
Dette problem kan angribes på flere forskellige måder.
i Ofte vælger man at transformere data så man har at gøre med et lineært problem. Antages f (K ) =
bK a finder vi ln f (K ) = ln b + a ln K. Derfor beskriver f (K ) sammenhængen mellem K og S netop
hvis
ln f (K ) = ln b + a ln K
beskriver den lineære sammenhæng mellem ln K og ln S. Vi kan da transformere data ved at tage
logaritmen til både kropsvægten og stofskiftet og plotte disse punkter. Dette er gjort i Figur 3. Den
rette linje, som bedst beskriver data har ligningen (der er afrundet til to decimaler)
ln f (K ) = −0.94 + 0.75 ln K.
4
Figur 3: Transformerede data (ln K, ln S) og den bedste rette linje
Det ses, at denne linje tilnærmer de transformerede punkter godt og derfor vil vi forvente en potenssammenhæng i data. Ordet “godt” kræver noget forklaring, som involverer modelcheck og som vi
udskyder til senere i uddannelsen. Transformeres tilbage kan vi dermed sige, at potensfunktionen
f (K ) = e−0.94+0.75 ln K = 0.39K0.75
bedst beskriver data.
ii “Back to basics”. Vi kan opskrive den tilsvarende sum af afvigelser Q( a, b) mellem f (K ) = bK a
og de ikke transformerede datapunkter. Denne sum af afvigelser betragtes som en funktion af de to
parametre og minimum for denne funktion bestemmes. Det er der to ting at sige til: For det første
kendes de matematiske metoder til at bestemme minimum for funktioner af to variable endnu ikke.
For det andet, selv når disse metoder kendes, er det en teknisk overvældende opgave at bestemme
minimum “i hånden” når blot datamængderne som her er bare lidt store. Funktionen er i dette
eksempel
6
Q( a, b) =
∑
S j − bK ja
2
,
j =1
hvor punkterne (K j , S j ) er angivet i tabelen oven for.
iii Bruge en mere avanceret “ikke lineær regression” på data. Basalt set er dette en implementering af
metoden “Back to basics”. Maple indeholder et udvalg af sådane rutiner. Med notationen oven fra
skrives
5
with(Statistics);
with(plots);
ker := [0.029, 1.1, 3.9, 7.0, 33, 48];
ser := [0.027, 0.4, 1.0, 2.2, 5.8, 6.0];
Fit(b*t^a, ker, ser, t);
Resultatet er, at a = 0.61 og b = 0.61, således, at den potensfunktion, der ligger tættest på punkterne f (K ) = 0.61 K0.61 .
iv Vi har altså to forskellige potensfunktioner, 0.61K0.761 og 0.39K0.75 , som begge kunne siges at
beskrive data bedst muligt. Som nævnt vil man ofte foretrække funktionen
f (K ) = 0.39K0.75 ,
fordi denne viser en lineær sammenhæng i det transformerede data. Det skal bemærkes, at parametrene a = 0.75 og b = 0.39 ikke minimerer funktionen
6
Q( a, b) =
∑
S j − bK ja
2
.
j =1
(Det gør derimod a = 0.61 og b = 0.61.) Derimod minimerer a = 0.75 og b = 0.39 funktionen
6
( a, b) →
∑
ln S j − ln b − a ln K j
2
.
j =1
v Ved at foretage en anden transformation af data fra Eksempel 2.2 kunne man undersøge en eventuel
eksponentiel sammenhæng i data: antages, at f (K ) = baK har vi ln f (K ) = ln b + K ln a hvorved
vi ville forvente at de transformerede punkter (K, ln S) ville tilnærmes godt af en ret linje. Det er
ikke tilfældet.
vi Inden for statistik er Maple iøvrigt ikke så udbredt som fx R-systemet. Det introduceres senere i
uddannelsen.
Bemærkning 2.1 Det er muligt at se lidt med, når Maple finder svarene til de opgaver man giver.
Dette gøres ved infolevel og skrives
infolevel[all] := 5;
Fit(b*t^a, ker, ser, t);
infolevel[all] := 0;
viser output, at Maple bruger NonlinearFit og LSSolve. Kort forklaret finder Maple først ud af, at der
skal anvendes en ikke linæer kurvetilpasning og dernæst metoden “Least Squares Solve”.
6
3
3.1
Modellering af hormesis
Hvad er hormesis?
Hormesis kan beskrives som et fænomen, hvor en lille dosis af “gift” eller anden form for biologisk stress til et biologisk system faktisk ikke er skadeligt for systemet men tværtimod kan være
gavnligt. Systemets respons på små doser er positiv, mens responsen på større doser er negativ.
Forskere er uenige om hvor udbredt og vigtigt fænomenet er. Fænomenet er formuleret relativt
bredt, som det fremgår nedenfor. Det turde også fremgå, at fænomenet i visse forklædninger
er ret kontroversielt.
• Indtagelse af små mængder alkohol kan have en positiv effekt på kredsløbet, mens indtagelse af større doser ikke har nogen gavnlig fysisk effekt!
• Strålingshormesis. Dette er en tese om, at små (vedvarende) doser af ioniserende stråling
er gavnlig for organismen. Den lave stråling kan stimulere organismens “reparationsmekanismer”, som derved kan beskytte mod forskellige former for sygdom. Laboratorieforsøg synes at bekræfte dette, men det har ikke været muligt at underbyggge teorien
med tilstrækkeligt data, og man opererer derfor generelt ikke med hormesis i forbindelse
med stråling.
Det kan nævnes, at byggematerialerne i forskellige boligblokke i Taiwan uoverlagt indeholdt den radioaktive isotop Cobolt-60. Undersøgelser viste imidlertid, at cancerdødsrater var meget mindre end i befolkningen generelt. Utroligt! og også temmeligt usandsynligt! Undersøgelserne tog nemlig ikke i betragtning, at beboerene i boligblokken
var relativt unge sammenlignet med Taiwans befolkning generelt. Faktisk viser senere
undersøgelser en signifikant vækst i fx leukæmi.
Neden for skal vi beskæftige os med hormesis inden for økotoksikologi. Pointen er ikke et
korstog for eller imod hormesis, men snarere at forsøge at beskrive forskellige modeller for
dosis-respons. Disse modeller indeholder parameter, der kan beskrive en hormetisk effekt, men
også kan udtale sig om hvorvidt en hormetisk effekt er signifikant. Vi giver en ramme for senere
at kunne undersøge, vha statistiske metoder, hvornår der synes at være noget i data, der tyder
på hormesis.
3.2
Økotoksikologi. Matematiske modeller
Inden for økotoksikologi undersøges fx planters respons når de udsættes for doser af forskellige toksiner. Man kan fx observere sammenhængen mellem koncentration af tilsat toksin og
bladlængde i en bestemt type plante.
For nogle kombinationer af planter og toksiner genkendes dosis-responskurverne som “logistiske” kurver af formen
d−c
y( x ) = c +
.
1 + eb(ln x−ln g)
7
Figur 4: En dosis-respons kurve
(Faktisk er disse kurver sammensætningen af den naturlige logaritme og en logistisk funktion
med bærekapacitet d − c og vækstrate −b.) Disse funktioner er aftagende og det illustrerer, at
større dosis medfører mindre respons. For andre kombinationer synes logistiske kurver ikke at
kunne beskrive sammenhængen mellem dosis og respons. Se bare på Figur 4, hvor responsen
af Cucumis sativus er studeret ved forskellige doser af toksinet parthenin.
Inden for Risk-management er man interesseret i niveauer, hvor p % af forskellen mellem kontrolresponsen (uden tilsætning af toksin) og den asymptotiske respons (tilsætning af meget
store doser toksin) er opnået. Disse niveauer kaldes EDp-niveauer og særligt ED50-niveauet
tjener som rettesnor for hvor meget toksin man kan nøjes med at give. I en kompleks situation, hvor flere toksiner og planter spiller sammen kan man måske komme ud for at et ED50
niveau for en plante bevirker en hormetisk effekt for andre planter og det kan være vanskeligt
at afgøre om det er optimalt at holde sig til ED50-niveauet. Bl.a. derfor er det vigtigt at kunne
afsløre hormetisk effekt.
En modificeret logistisk model, som tager højde for hormesis blev introduceret af Brain og
Cousens i 1989 (se [2]). Modellen er en responskurve af formen
y( x ) = c +
d−c+ fx
,
1 + eb(ln x−ln g)
som kun afviger fra en logistisk model ved inklusion af parameteren f .
Sætning 3.1 (Brain og Cousens model) Lad b, c, d, f og g være positive parametre og antag at b > 1
samt at d > c. Funktionen
d−c+ fx
, x>0
y( x ) = c +
1 + eb(ln x−ln g)
vokser fra værdien
lim y( x ) = d
x →0+
til et maksimum, hvorefter den aftager mod værdien
lim y( x ) = c.
x →∞
Set fra en statistisk vinkel kan denne model nogle gange være problematisk, idet et estimat
af parameteren f i visse tilfælde er negativt. Man kan vise, at den kurve der bedst beskriver
8
data da har en dal og ikke en top ved lave doser. Det er muligt at foretage parameterestimater
under antagelse af, at f er ikke-negativ, men det er ikke altid sket. Tilsvarende, hvis b er mindre
end 1 så beskriver kurven på ingen måde en dosis-respons kurve, og her vil estimater under
antagelsen b > 1 være nødvendige. Antagelser om både f og b kan lede til en mindre robust
model. Brain og Cousens model kan dog ofte bruges og den indgår i studier af hormesis i
blandinger, se fx [1]. Modellen er matematisk set relativt simpel, hvilket gør den nemmere at
bruge.
Cedergren, Ritz og Streibig (se [3]) har introduceret en mere robust, men også mere kompliceret
model for dosis-respons kurver. En af pointerne i dette skrift er at undersøge denne model.
Definition 3.1 (Cedergreen, Ritz og Streibig’s model) For positive parametre a, b og g defineres
−a
y( x ) = c +
−a
d − c + f e− x
d − c + f e− x
=
c
+
.
1 + ( x/g)b
1 + eb(ln x−ln g)
(De øvrige parametre er reelle.)
Forelagt data for en specifik plante og en specifik gift er det muligt statistisk at teste om parameteren f kunne være lig med 0. Lidt mere præcist testes det, hvor sandsynligt det ud fra data
er at f er forskellig fra 0. Man taler om, at parameteren f er signifikant (har en betydning).
Statistiske tests baseret på denne eller Brain og Cousens model et redskab for at afgøre om
hormesis optræder.
3.3
Egenskaber ved modelfunktionen. Parametrenes rolle
I dette afsnit undersøges hvordan de forskellige parametre påvirker responsfunktionen. Maple
bruges til at visualisere og hjælpe med (symbolske) beregninger. Funktionen er (som i Definition 3.1)
−a
−a
d − c + f e− x
d − c + f e− x
y( x ) = c +
=
c
+
,
1 + ( x/g)b
1 + eb(ln x−ln g)
hvor a, b og g antages positive. Vi noterer følgende simple egenskaber:
i Der gælder, at
lim y( x ) = d.
x →0
Parameteren d beskriver altså kontrolresponsen.
ii Tilsvarende repræsenterer parameteren c den asymptotiske respons (altså responsen ved
meget høje doser), idet der gælder
lim y( x ) = c.
x →∞
I den videre analyse vil vi antage, at der gælder d ≥ c, altså, at der for store doser er en toksisk
virkning. De øvrige parametre afspejler også egenskaber ved kurven, omend sammenhængen
mellem parametre og egenskaber ikke er helt simpel.
9
iii Antages, at f er positiv gælder der, at summen af parametrene d og f er en øvre grænse
−α
for den maksimale hormetiske effekt: Faktisk har vi e− x ≤ 1 og da eb(ln x−ln g) > 0 gælder
dermed
−a
y( x ) = c +
d − c + f e− x
d−c+ f
≤ c+
≤ c + d − c + f = d + f.
1 + eb(ln x−ln g)
1 + eb(ln x−ln g)
Det er vigtigt at pointere, at der er tale om en øvre grænse. Maksimum for funktionen y
kan i princippet udtrykkes vha alle de ingående parametre, men udregningerne bliver
hurtigt umulige: Den afledede kan udregnes vha Maple:
dr := x -> c+(d-c+f*exp(-1/x^a))/(1+exp(b*(ln(x)-ln(g))))
diff(dr(x),x)
Resultatet er
−a
y0 ( x ) =
x a +1
−a
f ae− x
(d − c + f e−x )beb(ln x−ln g)
−
.
2
1 + eb(ln x−ln g)
x 1 + eb(ln x−ln g)
iv Størrelsen af b har bl.a. betydning for hvor hastigt kurven aftager efter den maksimale
hormetiske effekt. Størrelsen af a har betydning for hvornår den maksimale hormetiske
effekt opnås. Se Opgave C.
Man kunne ledes til at tro, at funktionen i lighed med Brain og Cousens funktion vokser fra
værdien i x = 0 mod et maksimum, hvorefter den aftager mod grænseværdien for x gående
mod uendeligt. Dette er imidlertid ikke tilfældet, idet der gælder
lim y0 ( x ) x1−b = −(d − c)be−b ln g < 0.
x →0+
Heraf følger, at y( x ) er en aftagende funktion for x tæt på 0. Denne egenskab kunne synes at
være uheldig, men oftest er der slet ikke data ekstremt tæt på x = 0, og faktisk synes noget af
det data, der er tilgængeligt at “bekræfte” denne tendens.
Man kan endvidere vise, at y( x ) er aftagende for alle tilstrækkeligt store værdier af x. Disse
resultater udnytter størrelsesorden af eksponentialfunktioner og potensfunktioner; se Opgave
F.
3.4
Undersøgelse af EDp-niveauer
Et EDp-niveau er en dosis hvor p % af forskellen d − c mellem kontrolresponsen og den asymptotiske respons er opnået. Et sådant niveau x = EDp skal altså opfylde
−a
d − c + f e− x
p c+
=
c
+
1
−
.
100
1 + eb(ln x−ln g)
EDp-niveauer er interessante for koncentrationer, som er større end hvor den maksimale hormetiske effekt opnås. Til højre for maksimumspunktet for y er y en aftagende funktion og dette
sikrer, at der for hvert p er netop et tilhørende niveau EDp.
10
Figur 5: Dosis-respons med logaritmisk førsteakse. Parameterværdierne i Definition 3.1 er a =
0.25, b = 0.48, c = 0, d = 11, f = 416 og g = 0.014
Parameteren g er en nedre grænse for ED50 niveuaet, hvis f ≥ 0. Det betyder, at y( g) >
(c + d)/2, som indses på følgende måde:
−a
y( g) = c +
d − c + f e− g
d − c + f e− g
=
c
+
2
1 + eb(ln g−ln g)
−a
> c+
d−c
c+d
=
.
2
2
Til numerisk bestemmelse af EDp-niveauer anvendes generelle metoder for nulpunktsbestemmelse, som fx bisektionsmetoden, hvor et interval successivt midterdeles alt efter funktionens
(forskellige) fortegn i endepunkterne af intervallet.
3.5
Dosis-respons-grafik
Man ser, at x-aksen på den viste graf fra artiklen af Cedergreen, Ritz og Streibig er anderledes
end sædvanligt: aksen er logaritmisk. Se Figur 5. Data er observeret ved forskellige koncentrationer, hvor forholdet mellem en given koncentration og den næste er 10. Maple har indbygget
mulighed for at angive, at man ønsker en eller begge akser angivet i logaritmisk skala. Fx giver
følgende kode grafen for titalslogaritmen i Figur 6.
plot([log[10](x)], x=1..100, axis[1]=[mode=log, tickmarks=[1, 10, 100]])
3.6
Statistisk analyse af kurverne
I dette afsnit vil vi undersøge, hvordan Maple fungerer med kurvetilpasning og modellen af
Cedergren, Ritz og Streibig. Det er ikke altid, at Maple er lige pålidelig og det giver anledning
til implementering af en mere stabil og specialiseret rutine der kan fitte en generel klasse af
dosis-responskurver. En sådan rutine er udviklet af Ritz og Streibig, ikke til Maple men til R.
Se www.r-project.org under “Packages” og “drc”. Se endvidere www.bioassay.dk.
11
Figur 6: Et plot af titalslogaritmen med logaritmisk førsteakse
Dosis-responskurverne indeholder mange parametre. For at undgå overparametrisering foreslår Cedergreen, Ritz og Streibig, at parameteren a sættes til 0.25, 0.5 eller 1.0 alt efter hvonår en
mulig hormetisk effkt ser ud til at optræde.
4
Indlæsning af data i Maple
Her beskrives kort, hvilke (simple) muligheder Maple har i forbindelse med indlæsning af data
fra en fil. Der er også mere avancerede muligheder, der byder på større fleksibilitet mht format,
men dem vil vi forsøge at undgå.
• ImportMatrix. Filnavnet angives som en tekststreng ("filnavn") med angivelse af den
præcise placering af filen. Der er en række options, såsom overspringelse af linjer, angivelse af skilletegn og hvilken datatype der sal indlæses i. Eksempel:
ImportMatrix("/home/henrikp/undervisning/natit/byg.csv",
source = delimited, delimiter = ",", datatype = float[4], skiplines = 1)
• ImportData()– uden parametre. Her bliver man ført gennem en guide “Data Import Assistant”, hvorved man kan bestemme format osv.
• readdata Her indlæses data i søjler adskilt af blanktegn. Eksempel:
readdata("/home/henrikp/undervisning/natit/blanding.txt", [integer,
string, float, float])
12
5
Opgaver
A I filen byg.csv (findes på kursushjemmesiden) findes observerede data hvor en bygsort
Hordeum vulgare udsættes for herbicidet metsulfuron-methyl i forskellige doser. Responsen
måles i bladlængden i cm og herbicidet måles i µg/L.
(a) Undersøg indholdet i filen byg.csv fra kursushjemmesiden og download den.
(b) Indlæs data fra byg.csv i Maple fx vha
data := ImportMatrix("byg.csv", source = delimited,
delimiter = ",", datatype = float[4], skiplines = 1)
data er dermed en matrix (talskema) med 42 rækker og 2 søjler (skiplines bevirker, at der ses bort fra første linje i filen, der indeholder navnene på data). Brug evt
?ImportMatrix.
Check, at data er korrekt indlæst ved fx at skrive seq(data[k,1],k=1..42).
(c) Plot punkterne i et almindeligt xy-koordinatsystem (med plot(data, style=point)).
Kan man tydeligt se alle punkter?
(d) Bestem vha lineær regression den rette linje, der ligger tættest på punkterne. [Vink:
brug LinearFit.]
Bestem endvidere den logistiske funktion af formen
y ( x ) = a1 +
a2
,
1 + a 3 x a4
som bedst beskriver sammenhængen i data. [Vink: brug Fit]
Indtegn endeligt punkterne og de to regressionskurver i et koordinatsystem, hvor xaksen er logaritmisk. Kommentér denne graf. Hvor godt beskriver kurverne data?
(e) Plot datapunkterne igen som punkter i et nyt koordinatsystem, hvor x-aksen er logaritmisk.
Plot endvidere
−0.25
d − c + f e− x
y( x ) = c +
1 + eb(ln x−ln g)
i samme koordinatsystem for forskellige værdier af parametrene, og overbevis dig selv
om, at det ikke er lige til at finde de bedst mulige parametre. [Vink: brug eventuelt
numpoints=1000 i plottet for at få en mere glat kurve.]
(De optimale parameterværdier for a = 0.25 er b = 0.48, c = 0, d = 11, f = 416 og
g = 0.014.)
B Betragt funktionen
−a
d − c + f e− x
1 + eb(ln x−ln g)
med parameterværdierne a = 0.25, b = 0.48, c = 0, d = 11, f = 416 og g = 0.014. Lav
en funktionsundersøgelse af denne funktion vha Maple og bestem intervaller hvor funktionen aftager og intervaller hvor den vokser. [Vink: Det kan være nødvendigt med en meget
detaljeret graf tæt ved x = 0. Brug eventuelt også mode=log og fx numpoints=1000.]
y( x ) = c +
13
Bestem envidere ED50-niveauet for denne funktion. [Vink: fsolve.]
C I denne opgave undersøges, hvordan a og b indvirker på opførslen af grafen for y( x ) i Definition 3.1.
(a) Lad b = 2, c = 0, d = 1, g = 10 og f = 1. Undersøg vha Maple, hvad der ser ud til
at ske med den dosis, der resulterer i maksimal respons når a vokser. [Forslag: plot for
forskellige værdier af a i samme koordinatsystem; lav eventuelt en animation.]
(b) Lad a = 1, c = 0, d = 1, g = 10 og f = 1. Undersøg vha Maple, hvad der ser ud til at
ske med kurvens hældning efter maksimumsdosis med stigende b.
D Antag, at alle parametre er positive. Vis at y( x ) fra Definition 3.1 opfylder
c+
d−c
d−c+ f
≤ y( x ) ≤ c +
,
b
(
ln
x
−
ln
g
)
1+e
1 + eb(ln x−ln g)
og slut heraf, at kurven ligger mellem to logistiske kurver (som beskrevet i Afsnit 3.2).
E Lad nu a = 1, b = 2, c = 0, d = 1 og g = 1, men lad f være ukendt. Vis vha Maple, at
grafen for y( x ) i Definition 3.1 har en top for f = 10, er logistisk for f = 0 og har en dal for
f = −10.
F Eftervis, at y( x ) i Definition 3.1 er aftagende for alle x nær nul hvis der gælder d − c > 0 og
for alle tilstrækkeligt store x hvis der gælder d − c + f > 0.
G Bevis Sætning 3.1.
Litteratur
[1] R. G. Belz, N. Cedergreen and H. Sørensen, Hormesis in Mixtures - Can it be predicted?,
Science of the total Environment 404: 77–87, 2008.
[2] P. Brain and R. Cousens, An equation to describe dosis-responses where there is stimulation of
growth at low doses, Weed Res 29: 93–96, 1989.
[3] N. Cedergreen, C. Ritz and J.C. Streibig, Improved Empirical Models Describing Hormesis,
Environmental Toxicology and Chemistry 24: 3166–3172, 2005.
[4] T. V. Pedersen og H. L. Pedersen, Matematik og databehandling: Noter om matematik, Academic books, 2008.
14