data-raw/DATASET.R

In kyuhank/TMBworkshopPKNU:

## code to prepare `DATASET` dataset goes here

##################################################################################################################################
############################################ linreg data #########################################################################
##################################################################################################################################

linreg_pseudo_data=function(sample_size, alpha, beta, sd ) {

  set.seed(2)

  x=runif(sample_size, 1, 10)
  y=alpha+beta*x+rnorm(length(x), mean=0, sd=sd)

  return(list("x"=x, "y"=y))
}

sample_size=50
alpha=1
beta=2
sd=1.5

linreg_simul=linreg_pseudo_data(sample_size=sample_size, alpha=alpha, beta=beta, sd=sd)

linreg_pseudo=list("data"=linreg_simul, "input_pars"=list("alpha"=alpha, "beta"=beta, "sd"=sd))

usethis::use_data( linreg_pseudo, overwrite=T)

################################################################################################################################
############################################ SSPM data #########################################################################
################################################################################################################################

## Korean mackerel data
mackerel <- read.csv('mackerel_data.csv', header=F)
colnames(mackerel) <-c('yr','yield', 'CPUE')
usethis::use_data(mackerel)


## Simulated data
SSSP_simul_data2 <- function (seed,
                              data_size,
                              parameters_input=list("a"=a_true,
                                                    "b"=b_true,
                                                    "K"=K_true,
                                                    "r"=r_true,
                                                    "q"=q_true,
                                                    "sigma_pr"=sigma_pr_true,
                                                    "Exploit_min"=exp_lower,
                                                    "Exploit_max"=exp_upper)) {
  set.seed(seed)

  ndata=data_size

  ##bounded random walk for exploitation fractions
  x <- cumsum(sample(c(-1, 1), ndata, TRUE))
  E_true=((x-min(x)) /max(x - min(x) )* ( parameters_input[["Exploit_max"]] - parameters_input[["Exploit_min"]] ) )+parameters_input[["Exploit_min"]]

  Bt_true=c()
  logPt_true=c()
  Ct_true=c()
  It_obs=c()
  It_true=c()

  for (i in 2:ndata) {
    Bt_true[1]=parameters_input[["b"]]* parameters_input[["K"]] * rlnorm(1, meanlog=0, sdlog= parameters_input[["sigma_pr"]] )
    Ct_true[1]=Bt_true[1]*E_true[1]
    It_obs[1]=Bt_true[1]* parameters_input[["q"]]  * rlnorm(1, meanlog=0, sdlog= parameters_input[["a"]] * parameters_input[["sigma_pr"]] )
    Bt_true[i]=( Bt_true[i-1]+Bt_true[i-1] * parameters_input[["r"]] * ( 1 - (Bt_true[i-1] / parameters_input[["K"]] ) ) - Ct_true[i-1] )  * rlnorm(1, meanlog=0, sdlog= parameters_input[["sigma_pr"]] )
    Ct_true[i]=Bt_true[i]*E_true[i]
    It_obs[i]=Bt_true[i] * parameters_input[["q"]]  *rlnorm(1, meanlog=0, sdlog= parameters_input[["a"]] * parameters_input[["sigma_pr"]] )
  }

  It_true=Bt_true*parameters_input[["q"]]
  logPt_true=log(Bt_true/parameters_input[["K"]])

  names(Ct_true) <-paste("Y",1:data_size, sep='')
  names(Bt_true) <-paste("B",1:data_size, sep='')
  names(logPt_true) <-paste("logP",1:data_size, sep='')
  names(It_true) <-paste("I",1:data_size, sep='')
  names(It_obs) <-paste("I",1:data_size, sep='')

  return(list("Yt"=Ct_true,
              "It"=It_obs,
              "Bt"=Bt_true,
              #"true"=parameters_input
              #"It_true"=It_true,
              "logPt"=logPt_true)
         )
  }

seed=2

a=1
b=0.7
sigma_pr=0.1
K=50000
q=0.0001
r=0.3

parameters_simul=list("a"=a,
                      "b"=b,
                      "sigma_pr"=sigma_pr,
                      "K"=K,
                      "q"=q,
                      "r"=r,
                      "Exploit_min"=0.01,
                      "Exploit_max"=0.3 )

data_simul = SSSP_simul_data2 (seed=seed,data_size=43, parameters_input=parameters_simul)

Yt=data_simul$Yt[-43]
It=data_simul$It[-43]
logPt=data_simul$logPt
Bt=data_simul$Bt[-43]


sspm_pseudo=list("data"=list("Yt"=Yt, "It"=It),
                 "input_pars"=list("loga"=log(a), "logb"=log(b), "logsig_p"=log(sigma_pr), "logK"=log(K), "logq"=log(q), "logr"=log(r), "logPt"=logPt ),
                 "Bt"=Bt)

usethis::use_data(sspm_pseudo, overwrite=T)