# The n variable represents different values of 2N_{g}s
n <- 1

# The x variable is the parameter of the equation and represents u/v, the 
# mutation rate from "optimal" to "non-optimal" codon
x <- seq(0, 100, by=.0001)

y <- exp(n)/(exp(n)+x)

plot(x, y, type="l",
     xlab="Mutation rate toward non-optimal codon (u/v)",
     ylab="Equilibrium frequency of optimal codon (p̃)")

n <- 3
y <- exp(n)/(exp(n)+x)
lines(x, y=y, col="blue")

n <- 0
y <- exp(n)/(exp(n)+x)
lines(x, y=y, col="red")
# The n variable represents different values of 2N_{g}s
n <- 1

# The x variable is the parameter of the equation and represents u/v, the 
# mutation rate from "optimal" to "non-optimal" codon
x <- seq(0, 100, by=.0001)

y <- exp(n) / (exp(n)+x)

plot(x, y, log="x", type="l", 
     xaxt="n", 
     xlim=c(0.01, 100),
     xlab="Mutation rate toward non-optimal codon (u/v)",
     ylab="Equilibrium frequency of optimal codon (p̃)")

ticks <- seq(-2, 2, by=1)
labels <- sapply(ticks, function(i) as.expression(bquote(10^ .(i))))
axis(1, at=c(0.01, 0.1, 1, 10, 100), labels=labels)

n <- 3
y <- exp(n)/(exp(n)+x)
lines(x, y=y, col="blue")

n <- 0
y <- exp(n)/(exp(n)+x)
lines(x, y=y, col="red")
##--- Log-Log Plot  with  custom axes
lx <- seq(1, 5, length = 41)
yl <- expression(e^{-frac(1,2) * {log[10](x)}^2})
y <- exp(-.5*lx^2)

#op <- par(mfrow = c(2,1), mar = par("mar")-c(1,0,2,0), mgp = c(2, .7, 0))
plot(10^lx, y, log = "xy", type = "l", col = "purple",
     main = "Log-Log plot", ylab = yl, xlab = "x")
plot(10^lx, y, log = "xy", type = "o", pch = ".", col = "forestgreen",
     main = "Log-Log plot with custom axes", ylab = yl, xlab = "x",
     axes = FALSE, frame.plot = TRUE)
my.at <- 10^(1:5)
axis(1, at = my.at, labels = formatC(my.at, format = "fg"))
e.y <- -5:-1 ; at.y <- 10^e.y
axis(2, at = at.y, col.axis = "red", las = 1,
     labels = as.expression(lapply(e.y, function(E) bquote(10^.(E)))))
#par(op)


AlfonsoRReyes/RepDataPeerAssessment1 documentation built on May 5, 2019, 4:53 a.m.