Here is a similar example I'm using for my course, students seemed to enjoy
the simulation part at the end.

HTH,


Peter


*t<-100; alpha<-0.01; F<-0.6; N0<-80*

*r<-c(r0, numeric(t))*

*for (i in 1:t) *

*r[i+1]<- {*

*runif(1,0.2,1.8)*

*}*

*N<-c(N0, numeric(t))*

*for (i in 1:t) {*

*N[i+1]<-*

*(N[i]+r[i]*N[i]*(1-(N[i]*alpha))-F*N[i])*

*}*

*plot(0:t, N, type="l", xlab="time", ylab="fish population")*


*# Now repeat 100 times and some CI bands*


*t<-100*

*bevholt<-function(alpha=0.01, r0=1.2, F=0.6, N0=80, t=100) {*

*r<-c(r0, numeric(t))*

*for (i in 1:t) r[i+1]<- {*

*runif(1,0.4,1.4)*

*}*

*N<-c(N0, numeric(t))*

*for (i in 1:t) N[i+1]<- {*

*(N[i]+r[i]*N[i]*(1-(N[i]*alpha))-F*N[i])*

*}*

*return(N)*

*}*

*pop.sim<-bevholt()*

*plot(0:t, pop.sim, type="l", xlab="time", ylab="fish population")*


*sims<-replicate(10, {*

*output<-matrix(bevholt(), nrow=101)*

*},*

*simplify=TRUE)*

*matplot(0:t, sims, type="l")*


 *mean_sim<-apply(sims, 1, mean)*

*plot(0:t, mean_sim, type="l", ylim=c(15,80))*

*lower_ci<-mean_sim-apply(sims, 1, sd)*

*upper_ci<-mean_sim+apply(sims, 1, sd)*

*lines(lower_ci, lty='dashed')*

*lines(upper_ci, lty='dashed')*








On Mon, Feb 24, 2014 at 3:30 PM, Peter Solymos <solymos@ualberta.ca> wrote:

> Jeff,
>
> I am not sure why you need 100 random numbers for r and K, but if your goal
> is to get stochastic state-space model, you need to define the error term
> as a separate parameter and run the loop 100 times with the *same* fixed
> parameter values. When you do this, then you need to be aware of
> parametrization (e.g. Normal error on the log scale is what people usually
> use). Here is logistic Ricker state-space model:
>
> K <- 10000
> r <- 0.3
> sigma <- 0.2
> T <- 50
> N0 <- 10
>
> N <- numeric(T)
> N[1] <- N0
> for(i in 2:T) {
>     N[i] <- exp(r * (1 - N[i - 1]/K) + log(N[i - 1]) + rnorm(1, 0, sigma))
> }
> plot(N, type="b")
>
> Wrapped in a function, you get:
>
> Ricker <- function(r, K, sigma, T, N0) {
>     N <- numeric(T)
>     N[1] <- N0
>     for(i in 2:T) {
>         N[i] <- exp(r * (1 - N[i - 1]/K) + log(N[i - 1]) + rnorm(1, 0,
> sigma))
>     }
>     N
> }
>
> set.seed(1234)
> NN <- replicate(1000, Ricker(r=0.3, K=10^4, sigma=0.2, T=50, N0=10))
> matplot(NN, type="l", col=1, lty=1)
> lines(rowMeans(NN), col=2, lwd=2)
> abline(h=K, col=4, lty=2)
>
> HTH,
>
> Peter
>
>
> --
> Péter Sólymos
> 780-492-8534 | solymos@ualberta.ca | http://psolymos.github.com
> Alberta Biodiversity Monitoring Institute http://www.abmi.ca
> Boreal Avian Modelling Project http://www.borealbirds.ca
>
>
> On Sun, Feb 23, 2014 at 9:44 PM, Stratford, Jeffrey <
> jeffrey.stratford@wilkes.edu> wrote:
>
> > Hi everyone,
> >
> > I would like to create stochastic population models for an undergraduate
> > course.
> >
> > The goal would be to have students run models, record, results, change
> > parameters, and make inferences on changing the effects.
> >
> > I understand how to draw from distributions, I hit a knowledge wall with
> > loops.
> >
> > Here's what I have for a logistic model to predict N_t (not dN/dt)
> >
> > K <- rnorm(100, 10000, 400) # fluctuating environment
> > r <- rnorm(100, 1, 2) # variable reproductive potential
> > t <- 50
> > N0 <- 10 # start with 10 individuals
> >
> > for (i in 1:t) { N[,i+1] <- (N0*K)/(N0+((K-N0)^exp(-r*t)))}
> > plot(t, Nt)
> >
> > This doesn't work but I think I'm getting close.
> >
> > Any help would be much appreciated.
> >
> > Thanks,
> >
> > Jeff
> >
> > --
> > ********************************************************
> > Jeffrey A. Stratford, PhD
> > Department of Biology and Health Sciences
> > Wilkes University, PA 18766 USA
> > 570-332-2942
> > http://web.wilkes.edu/jeffrey.stratford/
> > ********************************************************
> >
> >         [[alternative HTML version deleted]]
> >
> > _______________________________________________
> > R-sig-ecology mailing list
> > R-sig-ecology@r-project.org
> > https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
> >
> >
>
>         [[alternative HTML version deleted]]
>
>
> _______________________________________________
> R-sig-ecology mailing list
> R-sig-ecology@r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
>
>


-- 

Peter Houk, PhD
Assistant Professor
University of Guam Marine Laboratory
http://www.guammarinelab.com/peterhouk.html
www.pacmares.com
www.micronesianfishing.com

	[[alternative HTML version deleted]]

