aridneptune aridneptune - 1 year ago 57
R Question

Plot density curve of mixture of two normal distribution

I am rather new to R and could use some basic help. I'd like to generate sums of two normal random variables (variance = 1 for each) as their means move apart and plot the results. The basic idea: if the means are sufficiently far apart, the distribution will be bimodal. Here's the code I'm trying:

x <- seq(-3, 3, length=500)
for(i in seq(0, 3, 0.25)) {
y <- dnorm(x, mean=0-i, sd=1)
z <- dnorm(x, mean=0+i, sd=1)
plot(x,y+z, type="l", xlim=c(-3,3))

Several questions:

  1. Are there better ways to do this?

  2. I'm only getting one PDF on my plot. How can I put multiple PDFs on the same plot?

Thank you in advance!

Answer Source

It is not difficult to do this using basic R features. We first define a function f to compute the density of mixture normal:

## `x` is an evaluation grid
## `dev` is deviation of mean from 0
f <- function (x, dev) {
  (dnorm(x, -dev) + dnorm(x, dev)) / 2

Then we use sapply to loop through various dev to get corresponding density:

## `dev` sequence to test
dev <- seq(0, 3, 0.25)
## evaluation grid; extending `c(-1, 1) * max(dev)` by 4 standard deviation
x <- seq(-max(dev) -4, max(dev) + 4, by = 0.1)
## density matrix
X <- sapply (dev, f, x = x)

Finally we use matplot for plotting:

matplot(x, X, type = "l", lty = 1)

enter image description here

More explanation

During sapply, x is not changed, while we pick up and try one element of dev each iteration. It is like

X <- matrix(0, nrow = length(x), ncol = length(dev))
for (i in 1:length(dev)) X[, i] <- f(x, dev[i])

matplot(x, X) will plot columns of X one by one, against x.