Mathematik für Biologiestudierende II¶
Sommersemester 2025
01.07.2025
© 2025 Prof. Dr. Rüdiger W. Braun
Klausur¶
- Sie können sich jetzt im Studierendenportal zur Klausur am 05.08.2025 anmelden.
- Die Klausur dauert 120 Minuten.
Termine¶
- 08.07.2025 letzte Vorlesung
- 15.07.2025 Probeklausur (im Vorlesungsformat)
- 05.08.2025 erste Klausur
- 16.09.2025 zweite Klausur
- Februar 2026 dritte Klausur
- möglicherweise wird Anfang September ein Ferientutoriumm angeboten: Schauen Sie auf die Webseite
Themen¶
- Lineare Modelle mit kategoriellen Daten
- ANOVA als lineares Modell
- Regression im exponentiellen Modell
- halblogarithmische Darstellung
In [1]:
import numpy as np
np.set_printoptions(legacy='1.21')
import pandas as pd
from scipy import stats
import seaborn as sns
sns.set_theme()
import statsmodels.formula.api as smf
Lineare Modelle mit kategoriellen Daten¶
Beispiel Ratten¶
Wir kommen zu dem Rattenbeispiel aus Lektion 23 zurück:
- kontaminiertes Gelände: fange 10 Ratten
- unbelastetes Vergleichsgelände: fange 10 Ratten
- für jede Ratte wird ihr Alter in Monaten und der Bleigehalt im Gewebe bestimmt
In [2]:
df = pd.read_csv('ratten.csv')
In [3]:
sns.lmplot(df, x='Alter', y='Belastung', hue='Gelände');
In [4]:
formel = 'Belastung ~ Alter + Gelände'
modell = smf.ols(formel, df)
res = modell.fit()
In [5]:
res.summary()
Out[5]:
| Dep. Variable: | Belastung | R-squared: | 0.673 |
|---|---|---|---|
| Model: | OLS | Adj. R-squared: | 0.635 |
| Method: | Least Squares | F-statistic: | 17.52 |
| Date: | Sun, 06 Jul 2025 | Prob (F-statistic): | 7.42e-05 |
| Time: | 17:01:42 | Log-Likelihood: | -63.935 |
| No. Observations: | 20 | AIC: | 133.9 |
| Df Residuals: | 17 | BIC: | 136.9 |
| Df Model: | 2 | ||
| Covariance Type: | nonrobust |
| coef | std err | t | P>|t| | [0.025 | 0.975] | |
|---|---|---|---|---|---|---|
| Intercept | 39.1728 | 5.166 | 7.583 | 0.000 | 28.273 | 50.072 |
| Gelände[T.unkontaminiert] | -11.0980 | 3.124 | -3.552 | 0.002 | -17.689 | -4.507 |
| Alter | 3.5490 | 0.617 | 5.752 | 0.000 | 2.247 | 4.851 |
| Omnibus: | 1.119 | Durbin-Watson: | 2.547 |
|---|---|---|---|
| Prob(Omnibus): | 0.572 | Jarque-Bera (JB): | 0.408 |
| Skew: | -0.346 | Prob(JB): | 0.815 |
| Kurtosis: | 3.102 | Cond. No. | 33.2 |
Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
Bestimmung des p-Werts¶
res.summarygibt den p-Wert für den Intercept aus- also für den Unterschied im Alter 0
- Trick: Verlegung des Nullpunkts
- Im Beispiel verlegen wird den Nullpunkt auf 8 Monate.
- Wir führen in der Tabelle also eine Spalte mit der Altersdifferenz zu 8 Monaten ein
In [6]:
df['Altersdifferenz'] = df.Alter - 8
df.head()
Out[6]:
| Alter | Belastung | Gelände | Altersdifferenz | |
|---|---|---|---|---|
| 0 | 10 | 63 | unkontaminiert | 2 |
| 1 | 12 | 67 | unkontaminiert | 4 |
| 2 | 6 | 55 | unkontaminiert | -2 |
| 3 | 6 | 42 | unkontaminiert | -2 |
| 4 | 11 | 73 | unkontaminiert | 3 |
In [7]:
formel2 = 'Belastung ~ Altersdifferenz + Gelände'
modell2 = smf.ols(formel2, df)
res2 = modell2.fit()
In [8]:
res2.summary()
Out[8]:
| Dep. Variable: | Belastung | R-squared: | 0.673 |
|---|---|---|---|
| Model: | OLS | Adj. R-squared: | 0.635 |
| Method: | Least Squares | F-statistic: | 17.52 |
| Date: | Sun, 06 Jul 2025 | Prob (F-statistic): | 7.42e-05 |
| Time: | 17:01:42 | Log-Likelihood: | -63.935 |
| No. Observations: | 20 | AIC: | 133.9 |
| Df Residuals: | 17 | BIC: | 136.9 |
| Df Model: | 2 | ||
| Covariance Type: | nonrobust |
| coef | std err | t | P>|t| | [0.025 | 0.975] | |
|---|---|---|---|---|---|---|
| Intercept | 67.5647 | 2.038 | 33.154 | 0.000 | 63.265 | 71.864 |
| Gelände[T.unkontaminiert] | -11.0980 | 3.124 | -3.552 | 0.002 | -17.689 | -4.507 |
| Altersdifferenz | 3.5490 | 0.617 | 5.752 | 0.000 | 2.247 | 4.851 |
| Omnibus: | 1.119 | Durbin-Watson: | 2.547 |
|---|---|---|---|
| Prob(Omnibus): | 0.572 | Jarque-Bera (JB): | 0.408 |
| Skew: | -0.346 | Prob(JB): | 0.815 |
| Kurtosis: | 3.102 | Cond. No. | 6.48 |
Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
- der p-Wert für den Einfluss des geländes bei acht Monate alten Ratten beträgt 0.002
- in der letzten Lektion hatten wir durch Betrachtung der Konfidenzintervalle der Prediction einen p-Wert von 0.02 gesehen
ANOVA als Lineares Modell¶
In [9]:
df = pd.read_csv("http://reh.math.uni-duesseldorf.de/~braun/bio2425/zitronen.csv")
In [10]:
df.head()
Out[10]:
| Vitamin_C_Gehalt | Land | |
|---|---|---|
| 0 | 494.5 | Spanien |
| 1 | 499.2 | Spanien |
| 2 | 494.3 | Spanien |
| 3 | 478.0 | Spanien |
| 4 | 500.1 | Spanien |
In [11]:
df.Land.value_counts()
Out[11]:
Land Spanien 8 Italien 8 Griechenland 8 Marokko 8 Indien 8 Name: count, dtype: int64
In [12]:
import statsmodels.stats.anova as smf_anova
In [13]:
formel = "Vitamin_C_Gehalt ~ Land"
modell = smf.ols(formel, df)
res = modell.fit()
In [14]:
smf_anova.anova_lm(res)
Out[14]:
| df | sum_sq | mean_sq | F | PR(>F) | |
|---|---|---|---|---|---|
| Land | 4.0 | 4378.44650 | 1094.611625 | 11.873758 | 0.000003 |
| Residual | 35.0 | 3226.56125 | 92.187464 | NaN | NaN |
- Der p-Wert ist derselbe wie in Lektion 17
- Die Nachkommastellen des p-Werts sehen wir wie folgt
In [15]:
tabelle = smf_anova.anova_lm(res)
In [16]:
tabelle['PR(>F)'].Land
Out[16]:
3.3733416696759452e-06
Die gemäß Bonferroni-Holm korrigierten p-Werte für die Paarvergleiche erhalten wir wie in Lektion 17
In [17]:
from statsmodels.sandbox.stats.multicomp import MultiComparison
In [18]:
muc = MultiComparison(df.Vitamin_C_Gehalt, df.Land)
In [19]:
res = muc.allpairtest(stats.ttest_ind, method='holm')
res[0]
Out[19]:
| group1 | group2 | stat | pval | pval_corr | reject |
|---|---|---|---|---|---|
| Griechenland | Indien | -4.9524 | 0.0002 | 0.0019 | True |
| Griechenland | Italien | 1.113 | 0.2845 | 0.5689 | False |
| Griechenland | Marokko | -3.5339 | 0.0033 | 0.0231 | True |
| Griechenland | Spanien | -1.9478 | 0.0718 | 0.2153 | False |
| Indien | Italien | 6.2008 | 0.0 | 0.0002 | True |
| Indien | Marokko | 0.3183 | 0.7549 | 0.7549 | False |
| Indien | Spanien | 3.3226 | 0.005 | 0.0302 | True |
| Italien | Marokko | -4.3312 | 0.0007 | 0.0055 | True |
| Italien | Spanien | -3.3042 | 0.0052 | 0.0302 | True |
| Marokko | Spanien | 2.2786 | 0.0389 | 0.1556 | False |
Das sind andere Zahlen als bei res.summary, weil res.summary die Treatments nur mit dem Default vergleicht und nicht untereinander
Regression im exponentiellen Modell¶
Beispiel Covid-Erkrankungen¶
In [20]:
df = pd.read_csv('corona.csv')
df.head()
Out[20]:
| Tag (im März) | Anzahl | |
|---|---|---|
| 0 | 3 | 38 |
| 1 | 4 | 52 |
| 2 | 5 | 109 |
| 3 | 6 | 185 |
| 4 | 7 | 150 |
In [21]:
ax = sns.scatterplot(data=df, x="Tag (im März)", y="Anzahl");
- Das Wachstum war exponentiell
- Es gab aber Schankungen durch unterschiedliche Verzögerungen der Berichte der Gesundheitsämter
Halblogarithmische Darstellung¶
Bei halblogarithmischer Darstellung
- ist die $x$-Achse linear skaliert: Gleiche absolute Zuwächse pro Längeneinheit
- ist die $y$-Achse logarithmisch skaliert: Gleiche relative Zuwächse pro Längeneinheit
- Das bedeutet: Der Logarithmus der Daten wird angezeigt, und die $y$-Achse wird entsprechend unterteilt
- Exponentiell wachsende Daten liegen bei halblogarithmischer Darstellung annäherend auf einer wachsenden Geraden, exponentiell fallende auf einer fallenden Geraden
In [22]:
ax.set_yscale('log')
ax.figure
Out[22]:
Exponentielles Modell vs. Lineare Regression¶
- Lineares Modell: in gleichen Zeitabständen gleiche absolute Zuwächse
- Exponentielles Modell: in gleichen Zeitabständen gleiche relative Zuwächse
- Biologische Wachstums- oder Abklingprozesse verlaufen meistens exponentiell
- Aufgabe der Regression im exponentiellen Modell ist es, bei Wachstumsprozessen die Verdoppelungszeit und bei Abklingprozessen die Halbwertszeit zu bestimmen
- Dies geschieht, indem man die Werte logarithmiert und dann deren lineare Regression berechnet
Regression im exponentiellen Modell¶
- $x$ die Zeit, $z$ Daten, die exponentiell wachsen (bzw. abklingen)
- Modellgleichung für Wachstumsprozess: $$ z = c \cdot e^{m\cdot x} $$
- logarithmierte Modellgleichung $$ y = \ln(z) = \ln(c) + m \cdot x $$
- bestimme diese Gerade durch lineare Regression
- wenn $m < 0$, dann Abklingprozess
In [23]:
df['logAnzahl'] = np.log(df.Anzahl)
df['Tag'] = df['Tag (im März)']
In [24]:
formel = 'logAnzahl ~ Tag'
modell = smf.ols(formel, df)
res = modell.fit()
In [25]:
res.summary()
Out[25]:
| Dep. Variable: | logAnzahl | R-squared: | 0.960 |
|---|---|---|---|
| Model: | OLS | Adj. R-squared: | 0.958 |
| Method: | Least Squares | F-statistic: | 458.1 |
| Date: | Sun, 06 Jul 2025 | Prob (F-statistic): | 9.25e-15 |
| Time: | 17:01:43 | Log-Likelihood: | -2.9636 |
| No. Observations: | 21 | AIC: | 9.927 |
| Df Residuals: | 19 | BIC: | 12.02 |
| Df Model: | 1 | ||
| Covariance Type: | nonrobust |
| coef | std err | t | P>|t| | [0.025 | 0.975] | |
|---|---|---|---|---|---|---|
| Intercept | 3.4410 | 0.151 | 22.728 | 0.000 | 3.124 | 3.758 |
| Tag | 0.2260 | 0.011 | 21.403 | 0.000 | 0.204 | 0.248 |
| Omnibus: | 0.968 | Durbin-Watson: | 1.597 |
|---|---|---|---|
| Prob(Omnibus): | 0.616 | Jarque-Bera (JB): | 0.816 |
| Skew: | -0.438 | Prob(JB): | 0.665 |
| Kurtosis: | 2.594 | Cond. No. | 34.1 |
Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
- m = 0.226
- b = 3.441
Die Regressionsgerade für die logarithmierten Daten ist $$ y = 0.226 \cdot x + 3.441 $$
In [26]:
tage = np.arange(3, 24)
gerade = 0.226*tage + 3.441
In [27]:
titel = "Die logarithmierten Daten zusammen mit ihrer Regressionsgerade"
ax2 = sns.scatterplot(x=df.Tag, y=df.logAnzahl)
sns.lineplot(x=tage, y=gerade)
ax2.set_title(titel);
In [28]:
titel = "Die exponentierte Regressionskurve zusammen mit den Ausgangsdaten in halblogarithmischer Darstellung"
sns.lineplot(x=tage, y=np.exp(gerade), ax=ax)
ax.set_title(titel)
ax.figure
Out[28]:
In [29]:
titel = "Die exponentierte Regressionskurve zusammen mit den Ausgangsdaten in linearer Darstellung"
ax.set_title(titel)
ax.set_yscale('linear')
ax.figure
Out[29]:
Halbwerts- bzw. Verdoppelungszeit¶
(aus Lektion 2)
- Modell eines Wachstumsprozesses $$ z = c \cdot e^{m\cdot x} $$
- Verdoppelungszeit $t$ bestimmt durch $$ e^{m\cdot t} = 2 $$
- Also $$ t = \frac{\ln 2}m $$
- Bei Abklingprozessen ist $m < 0$, dann ist $$ t = -\frac{\ln 2}m $$ die Halbwertszeit
Im Beispiel Covid
In [30]:
m = 0.226
In [31]:
t = np.log(2) / m
t
Out[31]:
3.067022922831616
Die Verdoppelungszeit betrug 3.07 Tage