output/script/self-assembly_pair.R

In Chang-Yu-Chang/MigrationCommunity: Microbial Consumer-Resource Model

library(tidyverse)
library(deSolve)
library(data.table)
library(micro.crm)

# Augument passed from command line
args = commandArgs(trailingOnly=TRUE)
treatment_name = args[1]
isolate1_freq = as.numeric(args[2])
isolate2_freq = as.numeric(args[3])

# Make regional species pool
parameter_load(parameter_text = paste0("parameter_script/", treatment_name, ".txt"))
pool_build(save_data = TRUE, data_name = "species_pool.Rdata")

# Make resource environment
resource <- setNames(c(1, rep(0, P-1)), paste0("R", sprintf("%03d", 1:P)))
assign("Rin", c(Rin1, rep(0, times = P-1)), envir = .GlobalEnv) # Supply for each resource. Set to 0

# Modify the `w` and `stoi` of fermenters.
w[1,1:(S_pool/function_group)] <- w_fermenter_P1
stoi[2,1,1:(S_pool/function_group)] <- stoi_fermenter_P1_P2

# Isolates ID
result_merged <- result_read(dir = "simulation_result/", pattern = treatment_name, threshold = threshold)

# Run
for (r in 1:20) {
  isolate_vector_r <- result_merged %>% filter(time == max(time), type == "X", replicate == r) %>% pull(variable)
  if (length(isolate_vector_r) == 1) next
  pair_df_r <- combn(isolate_vector_r,2) %>% t() %>% as.data.frame() %>% setNames(c("Isolate1", "Isolate2"))
  
  
  for (p in 1:nrow(pair_df_r)) {
    result = CR_model(c(resource, setNames(c(I*isolate1_freq, I*isolate2_freq), unlist(c(pair_df_r[p,])))),
             time_limit = time_limit, time_step = 1)
    # Write result file to txt
    parameter_write(result = result, 
                    parameter_list = parameter_list, 
                    file_name = paste0("simulation_result/pairs/", treatment_name, "_",
                                       sprintf("%03d", r), "_", isolate1_freq*100, "-", isolate2_freq*100, "_", 
                                       pair_df_r[p,1], "_", pair_df_r[p,2], ".txt"))
  }
  print (r)
}