Toadstools and Faeries
The purpose of this project is to bring together all of the piecesvectors, matrices, for loops, if statements, and functionsinto a single ecological modeling application.
Deep Background: each year toadstools in the forest produce faeries in the spring. In the fall, each faerie plants two toadstool seeds before dissolving in the final rays of the sun as it sets on the first day of winter. If there is a hard winter (because of climate change) then fewer of the toadstool seeds germinate in the spring. If it is a milder winter, more of the toadstools survive.
The purpose of this project is to model the 200year survival rate of toadstools given different probabilities of a hard winter. In order to do this, we will create the following functions.
 A function simYear to model a single year. The function returns the number of toadstools that survived the winter.
 A function runSim to model N years. The function returns a vector with the population at the end of earch year.
 A function runManySims to run Q simulations. The function returns a vector of the percent of nonzero (surviving) populations for each year of the simulation.
 A function calcProbSurvival that calculates the yearly survival probabilities given a list of different chances of a hard winter.
The final step will be to use the image function to view the resulting data. By adjusting some parameters of the toadstoolfaerie model, we can see the effects of changing different parameters.
Project Tasks

Model variables: place the following variable declarations at the top of your code file. These are the parameters for the model. By modifying their values, you can test the effects of different parameters on the model.
# model variables gbl_initial_pop_toadstools < 500 # initial population of toadstools gbl_chance_faerie < 0.8 # the chance a toadstool will create a faerie in the spring gbl_num_seeds_per_faerie < 2 # the number of seeds a faerie can carry gbl_chance_seed < 0.7 # the chance a seed will germinate gbl_hard_winter < 0.7 # the percent of seeds that survive a hard winter gbl_mild_winter < 0.9 # the percent of seeds that survive a mild winter gbl_max_toadstools < 1000 # the max number of toadstools (carrying capacity) gbl_min_toadstools < 10 # the minimum number of toadstools for a viable population # simulation variables gbl_NSims < 100 # number of simulations to execute gbl_NYears < 200 # number of years to simulation gbl_chance_of_hard_winter = c( 1/10.0, 1/9.0, 1/8.0, 1/7.0, 1/6.0, 1/5.0, 1/4.0, 1/3.0 )

Create the function simYear, which should have two parameters: num_toadstools, hard_winter_chance. It should calculate the number of toadstools remaining after the winter and return that value. The template for the function is below.
# calculates the number of toadstools after simulating a year # returns the population size at the end of the year simYear < function( num_toadstools, hard_winter_chance ) { # calculate the number of faeries created # code: assign to num_faeries the product of num_toadstools and gbl_chance_faerie # calculate the number of seeds that are planted # code: assign to num_seeds the product of the num_faeries and the gbl_num_seeds_per_faerie # figure out if it was a hard or mild winter and how many toadstools survive and germinate # code: assign to this_year a single value from runif # code: if this_year is less than hard_winter_chance # code: assign to spring_toadstools the product of num_seeds, gbl_hard_winter, and gbl_chance_seed # code: else # code: assign to spring_toadstools the product of num_seeds, gbl_mild_winter, and gbl_chance_seed # bounds check the number of toadstools # code: if spring_toadstools is greater than gbl_max_toadstools # code: assign to spring_toadstools the value in gbl_max_toadstools # code: if spring_toadstools is less than gbl_min_toadstools # code: assign to spring_toadstools the value 0 # code: return spring_toadstools }
Once you have written the above function, you can test it with the following code.
print( simYear( 500, 1.0 ) ) print( simYear( 500, 0.0 ) )
It should print out the values 392 and 504.

Create the function runSim, which should have two parameters: num_years, chance_hard_winter. The function should use the simYear function to simulate the specified number of years (num_years) and store each year's population size into a vector. The function should return a vector with num_years population values. The template is below.
# runs a simulation for the specified number of years and chance of a hard winter # returns a vector with the population size for each year of the simulation runSim < function(num_years, chance_hard_winter) { # code: assign to year_results a new double vector of size num_years # code: assign pop_toadstools the value in gbl_initial_pop_toadstools # code: for each year of the simulation (1:num_years) # code: assign to pop_toadstools the result of calling simYear with the arguments pop_toadstools and chance_hard_winter # code: assign to the next position in the year_results vector the value pop_toadstools (use your loop variable as the index) # code: return year_results }
Once you have written the above function, you can test it with the following code.
print( round(runSim( 10, 0.0 ),0) ) print( round(runSim( 10, 1.0 ),0) ) set.seed(42) print( round(runSim( 10, 0.2 ),0) )
You should get the following results
[1] 504 508 512 516 520 524 529 533 537 541 [1] 392 307 241 189 148 116 91 71 56 44 [1] 504 508 512 516 520 524 529 414 418 421

Create a function runManySims, which should have three parameters: num_sims, num_years, chance. The function should first execute num_sims simulations using the runSim function, storing the results in a matrix (one row for each simulation). It should then loop over each year and calculate the percentage of simulations with a population greater than 0, storing the result in a vector. The function should then return the vector of percentages.
# runs a group of simulations given number of simulations, number of years, and chance of a hard winter # returns a single vector with the probability of survival for each year runManySims < function( num_sims, num_years, chance ) { # calculate the matrix of num_sims v. num_years # code: assign to results a new matrix with num_years columns and num_sims rows # code: for i in the number of simulations (1:num_sims) code: assign to the ith row the result of calling runSim, passing in num_years and chance as the arguments # for each year, calculate the percent that is not zero # code: assign to prob_of_survival a new double vector of length num_years # code: for i in the number of years (1:num_years) # code: assign to the ith element of prob_of_survival the sum of the results in column that are greater than zero which is then divided by num_sims You can use the expression results[,i] > 0 to get a logical vector whose sum is the number of values in column i greater than 0. # code: return prob_of_survival }
Once you have written the above function, you can test it with the following code.
set.seed(42) print( round( runManySims( 100, 30, 0.5 ), 2 ) )
You should get the following results
[1] 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 [20] 0.99 0.99 0.99 0.99 0.98 0.96 0.94 0.87 0.81 0.73 0.64

Create a function calcProbSurvival, which should have three parameters: num_sims, num_years, chance_values. The function should execute the runManySims function for all of the values in the chance_values vector (chances of a hard winter), storing the results in a matrix. The function should return the matrix.
# runs a group of simulations, evaluating survival probabilities for each year # for a set of chance of hard winter values. # Returns a matrix with columns being chance of hard winter and rows being year calcProbSurvival < function( num_sims, num_years, chance_values ) { # code: assign to prob_survival a matrix with num_years rows and as many columns as the length of chance_values # code: assign to cur_col the value 1 # code: for each chance in chance_values # code: assign to colum cur_col of the prob_survival matrix the result of executing runManySims with the arguments num_sims, num_years, and chance # code: assign to cur_col the value of cur_col plus 1 # code: return prob_survival }
Once you have written the above function, run it using the following code.
all_results < calcProbSurvival( gbl_NSims, gbl_NYears, gbl_chance_of_hard_winter )
 Use the image function to plot the matrix all_results. Use 1:gbl_NYears as the x argument, use gbl_chance_of_hard_winter as the y argument, and use all_results as the z argument. Add appropriate x and y labels and a title.
Open R Studio or a text editor and Create a new R Script or R Markdown file (your choice) called project7. Put your name at the top of the file and save it to an appropriate working directory.
Report
Answer the following questions. Submit your answers as a plain text file or PDF in your handin directory. Put your name at the top of the file. No credit will be given for any other format or for a file without a name.
 If you change the the percent of toadstools that survive a hard winter (gbl_hard_winter) to 0.6, what happens?
 If you change the the percent of toadstools that survive a hard winter (gbl_hard_winter) to 0.8, what happens?
 Putting gbl_hard_winter back to 0.7, if you change the carrying capacity of the model (gbl_max_toadstools) to 200, what happens?
 Putting gbl_max_toadstools back to 1000, what happens if the chance of a toadstool making a faerie rises to 0.82?
Handin
Create a project7 folder inside the Private folder in your Courses directory. Put your R Script file, and your text/PDF file with your answers into the project7 folder. If you used R Markdown for the project, put your answers and the code all in the same file.