Effect_MOD_version2.R

getwd()

## [1] "C:/Users/hp/OneDrive/Desktop/Winter 2021/MOD3"

setwd("C:/Users/hp/OneDrive/Desktop/Winter 2021/MOD3")
library(drc)

## Loading required package: MASS

## 
## 'drc' has been loaded.

## Please cite R and 'drc' if used for a publication,

## for references type 'citation()' and 'citation('drc')'.

## 
## Attaching package: 'drc'

## The following objects are masked from 'package:stats':
## 
##     gaussian, getInitial

######################################
# ExerciseII; DAY1; EFFECT Modeling #
######################################
r <- 0.07
N <- 10000
tmax <- 72
EC <- 3.57
b <- -1.34
x <-rep(c(0), trunc(tmax))
x[30:35] <- 3
x[55:60] <- 5

res <- matrix(data=NA, nrow = tmax, ncol = 2)

for (t in 1:tmax){
  f <- 1/(1+((x[t]/EC)^b))
  dN <- (1-f)*r*N
  N <- N+dN
  
  res[t,1]<-t
  res[t,2]<-N
}

res

##       [,1]      [,2]
##  [1,]    1  10700.00
##  [2,]    2  11449.00
##  [3,]    3  12250.43
##  [4,]    4  13107.96
##  [5,]    5  14025.52
##  [6,]    6  15007.30
##  [7,]    7  16057.81
##  [8,]    8  17181.86
##  [9,]    9  18384.59
## [10,]   10  19671.51
## [11,]   11  21048.52
## [12,]   12  22521.92
## [13,]   13  24098.45
## [14,]   14  25785.34
## [15,]   15  27590.32
## [16,]   16  29521.64
## [17,]   17  31588.15
## [18,]   18  33799.32
## [19,]   19  36165.28
## [20,]   20  38696.84
## [21,]   21  41405.62
## [22,]   22  44304.02
## [23,]   23  47405.30
## [24,]   24  50723.67
## [25,]   25  54274.33
## [26,]   26  58073.53
## [27,]   27  62138.68
## [28,]   28  66488.38
## [29,]   29  71142.57
## [30,]   30  73921.46
## [31,]   31  76808.89
## [32,]   32  79809.11
## [33,]   33  82926.52
## [34,]   34  86165.70
## [35,]   35  89531.41
## [36,]   36  95798.60
## [37,]   37 102504.51
## [38,]   38 109679.82
## [39,]   39 117357.41
## [40,]   40 125572.43
## [41,]   41 134362.50
## [42,]   42 143767.87
## [43,]   43 153831.62
## [44,]   44 164599.84
## [45,]   45 176121.83
## [46,]   46 188450.35
## [47,]   47 201641.88
## [48,]   48 215756.81
## [49,]   49 230859.79
## [50,]   50 247019.97
## [51,]   51 264311.37
## [52,]   52 282813.17
## [53,]   53 302610.09
## [54,]   54 323792.80
## [55,]   55 332610.26
## [56,]   56 341667.84
## [57,]   57 350972.07
## [58,]   58 360529.67
## [59,]   59 370347.55
## [60,]   60 380432.78
## [61,]   61 407063.08
## [62,]   62 435557.49
## [63,]   63 466046.52
## [64,]   64 498669.77
## [65,]   65 533576.66
## [66,]   66 570927.03
## [67,]   67 610891.92
## [68,]   68 653654.35
## [69,]   69 699410.16
## [70,]   70 748368.87
## [71,]   71 800754.69
## [72,]   72 856807.52

par(mfrow = c(2, 1))
plot(res [,1], res [,2], type = "l", ylab = "N (cells/mL)", xlab = "Time (h)", main = "Population growth")
plot(res [,1], x, type = "l", ylab = "Conc", xlab = "Time (h)", main = "Chemical exposure")