Statistik mit R und RStudio - Ein Nachschlagewerk für Gesundheitsberufe

32 Diagramme plotten

Ich möchte Diagramme plotten!

32.1 Normalverteilung

Zum Plotten von Normalverteilungen kann man mit der Funktion plot() so vorgehen:

# erstelle Werte von -4 bis 4 in 0.005er-Schritten
x <- seq(-4, 4, by=0.005)

# plotte die Standardnormalverteilung
plot(x,dnorm(x))

# etwas hübscher
plot(x,dnorm(x), col="blue", type="l", xlab="x", ylab="f(x)", main="Standardnormalverteilungen")

# mit Mittelwert 2 und sd = 0,5
plot(x,dnorm(x,mean=2,s=0.5), col="darkblue", type="l", xlab="x", ylab="f(x)", main="Normalverteilungen")

# mit Mittelwert 2 und sd = 2
plot(x,dnorm(x,mean=2,s=2), col="darkorchid", type="l", xlab="x", ylab="f(x)", main="Normalverteilungen")

# erzeuge neue Werte von -3 bis 6
x <- seq(-3,6, by=0.005)

# Alles zusammen plotten
plot(x,dnorm(x,mean=2,s=0.5), col="blue", type="l", xlab="x", ylab="f(x)",main="Normalverteilungen")
lines(x,dnorm(x,mean=0,s=1), col="black")
lines(x,dnorm(x,mean=2,s=2), col="darkorchid")
text(0,.45,"N(0;1)")
text(2.8, 0.6, "N(2;0.5)", col="blue")
text(5, 0.1, "N(2;2)", col="darkorchid")

Mit ggplot kann es so aussehen:

# aktiviere ggplot
library(ggplot2)

# erzeuge Werte von -3 bis 3
x <- seq(-3,3, by=0.005)

# übergebe in ein data.frame 
df <- data.frame(x)

# ggplot erstellen
p <- ggplot(data=df, aes(x)) + 
                     xlim(-4,4) + ylim(0,0.5) + 
                     xlab("x") + ylab("Dichtefunktion") +
                     ggtitle("Normalverteilungen")
p + stat_function(fun=dnorm, args=(c(mean=0,sd=1)), colour="black")+ 
    annotate(geom="text", x=0, y=0.42, label="N(0;1)", color="black")

# erzeuge Werte von -3 bis 6
x <- seq(-3,6, by=0.005)

# übergebe in ein data.frame 
df <- data.frame(x)

# ggplot erstellen
p <- ggplot(data=df, aes(x)) + 
                     xlim(-3,6) + ylim(0, 0.8) + 
                     xlab("x") + ylab("Dichtefunktion") +
                     ggtitle("Normalverteilungen")
p + stat_function(fun=dnorm, args=(c(mean=0,sd=1)), colour="black")+ 
  annotate(geom="text", x=0, y=0.42, label="N(0;1)", color="black")+ 
  stat_function(fun=dnorm, args=(c(mean=2,sd=0.5)), colour="blue") +
  annotate(geom="text", x=3, y=0.6, label="N(2;0.5)", color="blue")+ 
  stat_function(fun=dnorm, args=(c(mean=2,sd=2)), colour="darkorchid") +
  annotate(geom="text", x=5, y=0.1, label="N(2;2)", color="darkorchid")

# erzeuge X-Werte
df    <- data.frame(x=seq(-3,3, by=0.005))
# berechne Y-Werte
df$y  <- dnorm(df$x)
# Setze Cuts
df$sd <- cut(df$x, breaks = c(-Inf, 1 *(-2:2), Inf))

# plotte als area
ggplot(df, aes(x, y, fill=sd)) + geom_area()

Erzeugen wir den Bereich von 1 Standardabweichung.

df <- data.frame(x=seq(-3,3, by=0.005))
df$y  <- dnorm(df$x)
df$sd <- cut(df$x, breaks = c(-1,1))
df <- rbind(df, data.frame(x=c(-1, 1), y=c(0,0),sd=c(-1,1)))
df$sd <- forcats::fct_explicit_na(df$sd, na_level="x") 

ggplot(df, aes(x, y, fill=sd)) + 
  ggtitle("Standardnormalverteilung", subtitle = "Bereich von 1 Standardabweichung") +
  geom_area() + theme(legend.position = "none") +
  geom_vline(xintercept=0, linetype="dotted")+
  geom_vline(xintercept=-1, linetype="dashed", col="blue", size=1)+
  geom_vline(xintercept=1, linetype="dashed", col="blue", size=1)+
  scale_fill_manual(name=c("x","(-1,1]"), values=c("skyblue", "snow3"))

Nun erzeugen wir den Bereich von 2 Standardabweichungen.

df    <- data.frame(x=seq(-3,3, by=0.005))
df$y  <- dnorm(df$x)
df$sd <- cut(df$x, breaks = c(-2,2))
df    <- rbind(df, data.frame(x=c(-2, 2), y=c(0,0),sd=c(-2,2)))
df$sd <- forcats::fct_explicit_na(df$sd, na_level="x") 

ggplot(df, aes(x, y, fill=sd)) + 
  ggtitle("Standardnormalverteilung", subtitle = "Bereich von 2 Standardabweichung") +
  geom_area() + theme(legend.position = "none") +
  geom_vline(xintercept=0, linetype="dotted")+  
  geom_vline(xintercept=-2, linetype="dashed", col="purple3", size=1)+
  geom_vline(xintercept=2, linetype="dashed", col="purple3", size=1)+
  scale_fill_manual(name=c("x","(-2,2]"), values=c("purple3", "snow3"))

Erzeugen wir die Fläche, in der 95% der Werte liegen, also von 1,96 Standardabweichungen.

df    <- data.frame(x=seq(-3,3, by=0.005))
df$y  <- dnorm(df$x)
df$sd <- cut(df$x, breaks = c(-1.96, 1.96))
df    <- rbind(df, data.frame(x=c(-1.96, 1.96), y=c(0,0),sd=c(-1.96,1.96)))
df$sd <- forcats::fct_explicit_na(df$sd, na_level="x") 

ggplot(df, aes(x, y, fill=sd)) + 
  ggtitle("Standardnormalverteilung", subtitle = "95% der Werte bei 1,96 sd") +
  geom_area() + theme(legend.position = "none") +
  geom_vline(xintercept=0, linetype="dotted")+
  geom_vline(xintercept=-2, linetype="dotted", col="purple3", size=0.5)+
  geom_vline(xintercept=2, linetype="dotted", col="purple3", size=0.5)+
  geom_vline(xintercept=-1.96, col="chartreuse4", size=0.5)+
  geom_vline(xintercept=1.96, col="chartreuse4", size=0.5)+
  scale_fill_manual(name=c("x","(-1.96,1.96]"), values=c("seagreen3", "snow3"))

Dieser Code kann für die Erzeugung der Grafiken von Tabelle Tabelle 36.1 und Tabelle 36.2 verwendet werden.

df    <- data.frame(x=seq(-3,3, by=0.005))
df$y  <- dnorm(df$x)
df$sd <- cut(df$x, breaks = c(-3, 1.4))
df    <- rbind(df, data.frame(x=c(-3, 1.4), y=c(0,0),sd=c(-3,1.4)))
df$sd <- forcats::fct_explicit_na(df$sd, na_level="x") 

ggplot(df, aes(x, y, fill=sd)) + 
  geom_area() + theme(legend.position = "none") +
  geom_vline(xintercept=0, linetype="dotted")+
  geom_vline(xintercept=1.4, col="blue", size=0.5)+
  ylab("") + scale_fill_manual(values=c("skyblue", "snow3"))

df    <- data.frame(x=seq(-3,3, by=0.005))
df$y  <- dnorm(df$x)
df$sd <- cut(df$x, breaks = c(-3, 1.4))
df    <- rbind(df, data.frame(x=c(-3, 1.4), y=c(0,0),sd=c(-3,1.4)))
df$sd <- forcats::fct_explicit_na(df$sd, na_level="x") 

ggplot(df, aes(x, y, fill=sd)) + 
  geom_area() + theme(legend.position = "none") +
  geom_vline(xintercept=0, linetype="dotted")+
  geom_vline(xintercept=1.4, col="blue", size=0.5)+
  ylab("") + scale_fill_manual(values=c("snow3", "skyblue"))

32.2 t-Verteilung

Die t-Verteilung kann mit ggplot geplottet werden.

# Erzeuge x-werte
df <- data.frame(x=seq(-3,3, by=0.005))

# Grundlegene Plotangaben
p <- ggplot(data=df, aes(x)) +
  # begrenze die Achsen
  xlim(-3,3) + ylim(0, 0.4) +
  # Achsen-Titel
  xlab("x") + ylab("Dichtefunktion") +
  # Plot-Titel
  ggtitle("t-Verteilungen", subtitle = "nach Freiheitsgraden")

  # t-Verteilung plotten
p + 
  stat_function(fun=dt, args=list(df=1), col="black") +
  # Textfeld hinzufügen
  annotate(geom="text", x=0, y=0.25, label="df=1", color="black")

Dem Plott können weitere Freiheitsgrade hinzugefügt werden.

# t-Verteilungen plotten
p + 
  stat_function(fun=dt, args=list(df=1), col="black") +
  annotate(geom="text", x=0, y=0.25, label="df=1", color="black")+
  stat_function(fun=dt, args=list(df=2), col="blue") +
  annotate(geom="text", x=0, y=0.37, label="df=2", color="blue")

Wenn die t-Werte bereits berechnet wurden, kann eine alternative Vorgehensweise so aussehen:

# berechne t-Werte für Freiheitsgrade 1 bis 5
x=seq(-3,3, by=0.005)
df <- data.frame(
    x,
    df1 = dt(x,df=1), 
    df2 = dt(x,df=2),
    df3 = dt(x,df=3),
    df4 = dt(x,df=4),
    df5 = dt(x,df=5)
)

# wandle ins Format "long table" um
df <- pivot_longer(df, cols=c(df1, df2, df3, df4, df5))

# grundlegende Ploteinstellungen
p <- ggplot(data=df, aes(x,value)) +
  xlim(-3,3) + ylim(0, 0.4) +
  xlab("x") + ylab("Dichtefunktion") +
  ggtitle("t-Verteilungen", subtitle = "nach Freiheitsgraden")

# plotte t-Verteilungen
p + geom_line(aes(col=name))+
  labs(col="Freiheitsgrade")

32.3 \(\chi^2\)-Verteilung

Die \(\chi^2\)-Verteilung kann mit ggplot geplottet werden.

# Erzeuge x-werte
x=seq(0,25, by=0.005)
df <- data.frame(  x,
  df01 = dchisq(x, df=1), 
  df05 = dchisq(x, df=5),
  df10 = dchisq(x, df=10),
  df15 = dchisq(x, df=15)
)

# erzeuge long-table
df <- pivot_longer(df, cols=c(df01, df05, df10, df15))

p <- ggplot(data=df, aes(x, value, fill=name)) +
  xlim(0,25) + ylim(0, 0.2) +
  xlab("x") + ylab("Dichtefunktion") +
  ggtitle("Chi^2-Verteilungen", subtitle = "nach Freiheitsgraden")
p + geom_line(aes(col=name,linetype=name))+
labs(col="Freiheitsgrade",linetype="")

Die Fläche unterhalb der Kurve kann mit geom_area() erzeugt werden. Für df=1 lautet der Aufruf:

# erzeuge Dummy-Werte
x=seq(0,25, by=0.005)
# überführe in Datenframe
df <- data.frame(  x, df01 = dchisq(x, df=1))
ggplot(data=df, aes(x, df01)) +
  xlim(0,25)  + coord_cartesian(ylim=c(0, 0.2)) +
  xlab("Chi^2-Wert") + ylab("Dichtefunktion") +
  ggtitle("Chi^2-Verteilung", subtitle = "mit 1 Freiheitsgrad") +   
  geom_line(col="#F8766D") +geom_area(fill="skyblue")

Ähnlich wie bei der Normalverteilung können wir die Fläche für bestimmte \(\chi^2\)-Werte einfärben. Bei einem Freiheitsgrad und \(\alpha=0,05\) ergibt sich ein kritischer \(\chi^2\)-Wert von 3,84. Der Flächenanteil unterhalb dieses Wertes lässt sich wie folgt darstellen:

# Dummy-Werte erzeugen
df <- data.frame(x=seq(0.005,25, by=0.005))
# Chi^2-Werte für df=1 erstellen
df$y  <- dchisq(df$x, df=1)
# Grenze bei 3.84 einziehen
df$sd <- cut(df$x, breaks = c(0.00, 3.84))
ggplot(df, aes(x, y, fill=sd)) + 
  geom_area() + theme(legend.position = "none") + 
  geom_line(aes(col="#F8766D")) +
  xlim(0,25) + coord_cartesian(ylim=c(0, 0.2)) +
  xlab("Chi^2-Wert") + ylab("Dichtefunktion") +
  ggtitle("Chi^2-Verteilung", subtitle = "mit 1 Freiheitsgrad")+
  geom_vline(xintercept=3.84, col="blue", size=0.5, linetype="dashed")+
  ylab("") + scale_fill_manual(values=c("skyblue", "snow3")) +
  annotate(geom="text", x=2, y=0.028, size=8,label="95%", color="blue")

32.4 Poisson-Verteilung

Die Poisson-Verteilung kann mit ggplot geplottet werden.

# Erzeuge x-Werte
x=seq(0,25)
df <- data.frame(  x,
                   l1 = dpois(x, 1), 
                   l2 = dpois(x, 2),
                   l4 = dpois(x, 4),
                   l9 = dpois(x, 9)
)

# erzeuge eine long-table
df <- pivot_longer(df, cols=c(l1, l2, l4, l9))

# plot vorbereiten
p <- ggplot(data=df, aes(x, value, fill=name)) +
  xlim(0,17) + ylim(0, 0.4) +
  xlab("x") + ylab("Dichtefunktion") +
  ggtitle("Poisson-Verteilungen", subtitle = "nach Lambda")
p + geom_line(aes(col=name,linetype=name))+
  labs(col="Lambda",linetype="")

## Warning: Removed 32 row(s) containing missing values (geom_path).

Da wir die Plot-Grundlagen in p gespeichert haben, können wir ergänzen:

p + geom_bar(aes(col=name,linetype=name), stat="identity", position="dodge")+
   geom_line(aes(col=name,linetype=name))+
  labs(col="Lambda",linetype="")

## Warning: Removed 36 rows containing missing values (geom_bar).

## Warning: Removed 32 row(s) containing missing values (geom_path).

Möchte man einen Polygonzug, müssen als Ausgangspunkt die Koordinaten (0,0) festegelegt werden.

x=seq(0,25)#, by=0.005)
df <- data.frame(  x=c(0,x),
                   l1 = c(0, dpois(x, 1)), 
                   l2 = c(0, dpois(x, 2)),
                   l4 = c(0, dpois(x, 4)),
                   l9 = c(0, dpois(x, 9))
)
df <- pivot_longer(df, cols=c(l1, l2, l4, l9))

p <- ggplot(data=df, aes(x, value, fill=name)) +
  xlim(0,17) + ylim(0, 0.4) +
  xlab("x") + ylab("Dichtefunktion") +
  ggtitle("Poisson-Verteilungen", subtitle = "nach Lambda")
p + geom_polygon(aes(col=name,linetype=name))+
  geom_line(aes(col=name),linetype="dashed")+
  labs(col="Lambda", fill="Lambda", linetype="Lambda")

## Warning: Removed 32 row(s) containing missing values (geom_path).