Mercurial > repos > greg > insect_phenology_model

comparison insect_phenology_model.R @ 47:d6482112bf35 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author greg
date Mon, 13 Nov 2017 10:53:31 -0500
parents d7e406de8094
children 66712cdec7c0
comparison
equal deleted inserted replaced
46:d7e406de8094 47:d6482112bf35
58 dnp <- temperature_data[row, 5] 58 dnp <- temperature_data[row, 5]
59 # Maximum temperature for current row. 59 # Maximum temperature for current row.
60 dxp <- temperature_data[row, 6] 60 dxp <- temperature_data[row, 6]
61 # Mean temperature for current row. 61 # Mean temperature for current row.
62 dmean <- 0.5 * (dnp + dxp) 62 dmean <- 0.5 * (dnp + dxp)
63 dd <- 0 # initialize degree day accumulation 63 # Initialize degree day accumulation
64 dd <- 0
64 if (dxp < threshold) { 65 if (dxp < threshold) {
65 dd <- 0 66 dd <- 0
66 } 67 }
67 else { 68 else {
68 # Initialize hourly temperature. 69 # Initialize hourly temperature.
237 # Photoperiod in the day. 238 # Photoperiod in the day.
238 photoperiod <- daylight_length_vector[row] 239 photoperiod <- daylight_length_vector[row]
239 temp.profile <- get_temperature_at_hour(latitude, temperature_data, daylight_length_vector, row, opt$num_days) 240 temp.profile <- get_temperature_at_hour(latitude, temperature_data, daylight_length_vector, row, opt$num_days)
240 mean.temp <- temp.profile[1] 241 mean.temp <- temp.profile[1]
241 dd.temp <- temp.profile[2] 242 dd.temp <- temp.profile[2]
242 dd.day[day] <- dd.temp 243 dd.day[row] <- dd.temp
243 # Trash bin for death. 244 # Trash bin for death.
244 death.vec <- NULL 245 death.vec <- NULL
245 # Newborn. 246 # Newborn.
246 birth.vec <- NULL 247 birth.vec <- NULL
247 248
266 else if (vec.ind[2] == 1 | vec.ind[2] == 2) { 267 else if (vec.ind[2] == 1 | vec.ind[2] == 2) {
267 death.prob = opt$nymph_mort * mortality.nymph(mean.temp) 268 death.prob = opt$nymph_mort * mortality.nymph(mean.temp)
268 } 269 }
269 else if (vec.ind[2] == 3 | vec.ind[2] == 4 | vec.ind[2] == 5) { 270 else if (vec.ind[2] == 3 | vec.ind[2] == 4 | vec.ind[2] == 5) {
270 # For adult. 271 # For adult.
271 if (day < day.kill) { 272 if (doy < day.kill) {
272 death.prob = opt$adult_mort * mortality.adult(mean.temp) 273 death.prob = opt$adult_mort * mortality.adult(mean.temp)
273 } 274 }
274 else { 275 else {
275 # Increase adult mortality after fall equinox. 276 # Increase adult mortality after fall equinox.
276 death.prob = opt$adult_mort * post.mort * mortality.adult(mean.temp) 277 death.prob = opt$adult_mort * post.mort * mortality.adult(mean.temp)
284 else { 285 else {
285 # Aggregrate index of dead bug. 286 # Aggregrate index of dead bug.
286 # Event 1 end of diapause. 287 # Event 1 end of diapause.
287 if (vec.ind[1] == 0 && vec.ind[2] == 3) { 288 if (vec.ind[1] == 0 && vec.ind[2] == 3) {
288 # Overwintering adult (previttelogenic). 289 # Overwintering adult (previttelogenic).
289 if (photoperiod > opt$photoperiod && vec.ind[3] > 68 && day < 180) { 290 if (photoperiod > opt$photoperiod && vec.ind[3] > 68 && doy < 180) {
290 # Add 68C to become fully reproductively matured. 291 # Add 68C to become fully reproductively matured.
291 # Transfer to vittelogenic. 292 # Transfer to vittelogenic.
292 vec.ind <- c(0, 4, 0, 0, 0) 293 vec.ind <- c(0, 4, 0, 0, 0)
293 vec.mat[i,] <- vec.ind 294 vec.mat[i,] <- vec.ind
294 } 295 }
350 birth.vec <- rbind(birth.vec, new.vec) 351 birth.vec <- rbind(birth.vec, new.vec)
351 } 352 }
352 } 353 }
353 354
354 # Event 2 oviposition -- for gen 1. 355 # Event 2 oviposition -- for gen 1.
355 if (vec.ind[2] == 4 && vec.ind[1] == 1 && mean.temp > 12.5 && day < 222) { 356 if (vec.ind[2] == 4 && vec.ind[1] == 1 && mean.temp > 12.5 && doy < 222) {
356 # Vittelogenic stage, 1st generation 357 # Vittelogenic stage, 1st generation
357 if (vec.ind[4] == 0) { 358 if (vec.ind[4] == 0) {
358 # Just turned in vittelogenic stage. 359 # Just turned in vittelogenic stage.
359 n.birth=round(runif(1, 2 + opt$min_clutch_size, 8 + opt$max_clutch_size)) 360 n.birth=round(runif(1, 2 + opt$min_clutch_size, 8 + opt$max_clutch_size))
360 } 361 }
411 if (vec.ind[3] >= (250 + opt$old_nymph_accum)) { 412 if (vec.ind[3] >= (250 + opt$old_nymph_accum)) {
412 # From young to old nymph, dd requirement met. 413 # From young to old nymph, dd requirement met.
413 current.gen <- vec.ind[1] 414 current.gen <- vec.ind[1]
414 # Transfer to old nym stage. 415 # Transfer to old nym stage.
415 vec.ind <- c(current.gen, 2, 0, 0, 0) 416 vec.ind <- c(current.gen, 2, 0, 0, 0)
416 if (photoperiod < opt$photoperiod && day > 180) { 417 if (photoperiod < opt$photoperiod && doy > 180) {
417 vec.ind[5] <- 1 418 vec.ind[5] <- 1
418 } # Prepare for diapausing. 419 } # Prepare for diapausing.
419 } 420 }
420 else { 421 else {
421 # Add 1 day in current stage. 422 # Add 1 day in current stage.
489 gen1 <- sum(vec.mat[,1] == 1) 490 gen1 <- sum(vec.mat[,1] == 1)
490 # Second generation. 491 # Second generation.
491 gen2 <- sum(vec.mat[,1] == 2) 492 gen2 <- sum(vec.mat[,1] == 2)
492 # Sum of all adults. 493 # Sum of all adults.
493 n.adult <- sum(vec.mat[,2] == 3) + sum(vec.mat[,2] == 4) + sum(vec.mat[,2] == 5) 494 n.adult <- sum(vec.mat[,2] == 3) + sum(vec.mat[,2] == 4) + sum(vec.mat[,2] == 5)
494 # Gen eration 0 pop size. 495
495 gen0.pop[day] <- gen0 496 # Generation 0 pop size.
496 gen1.pop[day] <- gen1 497 gen0.pop[row] <- gen0
497 gen2.pop[day] <- gen2 498 gen1.pop[row] <- gen1
498 S0[day] <- s0 499 gen2.pop[row] <- gen2
499 S1[day] <- s1 500
500 S2[day] <- s2 501 S0[row] <- s0
501 S3[day] <- s3 502 S1[row] <- s1
502 S4[day] <- s4 503 S2[row] <- s2
503 S5[day] <- s5 504 S3[row] <- s3
504 g0.adult[day] <- sum(vec.mat[,1] == 0) 505 S4[row] <- s4
505 g1.adult[day] <- sum((vec.mat[,1] == 1 & vec.mat[,2] == 3) | (vec.mat[,1] == 1 & vec.mat[,2] == 4) | (vec.mat[,1] == 1 & vec.mat[,2] == 5)) 506 S5[row] <- s5
506 g2.adult[day] <- sum((vec.mat[,1]== 2 & vec.mat[,2] == 3) | (vec.mat[,1] == 2 & vec.mat[,2] == 4) | (vec.mat[,1] == 2 & vec.mat[,2] == 5)) 507
507 508 g0.adult[row] <- sum(vec.mat[,1] == 0)
508 N.newborn[day] <- n.newborn 509 g1.adult[row] <- sum((vec.mat[,1] == 1 & vec.mat[,2] == 3) | (vec.mat[,1] == 1 & vec.mat[,2] == 4) | (vec.mat[,1] == 1 & vec.mat[,2] == 5))
509 N.death[day] <- n.death 510 g2.adult[row] <- sum((vec.mat[,1]== 2 & vec.mat[,2] == 3) | (vec.mat[,1] == 2 & vec.mat[,2] == 4) | (vec.mat[,1] == 2 & vec.mat[,2] == 5))
510 N.adult[day] <- n.adult 511
512 N.newborn[row] <- n.newborn
513 N.death[row] <- n.death
514 N.adult[row] <- n.adult
511 } # end of days specified in the input temperature data 515 } # end of days specified in the input temperature data
512 516
513 dd.cum <- cumsum(dd.day) 517 dd.cum <- cumsum(dd.day)
518
514 # Collect all the outputs. 519 # Collect all the outputs.
515 S0.rep[,N.rep] <- S0 520 S0.rep[,N.rep] <- S0
516 S1.rep[,N.rep] <- S1 521 S1.rep[,N.rep] <- S1
517 S2.rep[,N.rep] <- S2 522 S2.rep[,N.rep] <- S2
518 S3.rep[,N.rep] <- S3 523 S3.rep[,N.rep] <- S3