load("Data/All_simulations_combined/savedsims4.Rdata")
length(savedsims)
## [1] 128

The savedsims data contains ALL the data. It lists the simulation data output from QPress for each of the 128 model versions, where each element on the list is the output from one model version.

#the nodes in the models
nodes<-c("BigFish", "Bycatch", "Catch", "DeepPredators", "DeepPrey", "Depredation", "Downwelling", "Effort",   "Landings", "LiceMBPs", "LocalDepletion", "MoveOn", "PreRecruits", "Quota", "ShallowPrey", "SmallFish", "SST", "TempDepth")

We need to get the edge strength and node response data from across the different model versions listed in savedsims.

Edge strength (weight) data is found in the lists w and node responses in lists A for each model. In each of these lists there is a matrix containing the data for each of the 10000 simulations. Each row in the A matrices corresponds to one press perturbation, with the SST perturbation in row 17. The columns correspond to the response of each node to these press perturbations, in alphabetical order.

rPRec<-vector()
rSP<-vector()
rSF<-vector()
rDP<-vector()
rBF<-vector()
rC<-vector()
rL<-vector()
rLice<-vector() 
mv<-vector()
TtoSP<-vector()
sTSP<-vector()
sTPR<-vector()


# SST Press responses
for (j in 1:length(savedsims)){
  responsePRec<-0
  responseSP<-0
  responseSF<-0
  responseDP<-0
  responseBF<-0
  responseC<-0
  responseL<-0
  responseLice<-0
  TemptoSP<-0
  modver<-0
  strengthTSP<-0
  strengthTPR<-0

  
  
    for (i in 1:10000){
      
      #get the edge weights for the two outbound edges from the SST node
        strengthTSP[i]<-savedsims[[j]]$w[i, grep("^SST .* ShallowPrey$", colnames(savedsims[[j]]$w))] #col#14
        strengthTPR[i]<-savedsims[[j]]$w[i, grep("^SST .* PreRecruits$", colnames(savedsims[[j]]$w))] #SST-*PreRecruits, col#13
        
      #get the responses to the SST perturbation
      responseL[i]<-savedsims[[j]]$A[[i]][9, 17]
      responseBF[i]<-savedsims[[j]]$A[[i]][1,17]
      responseC[i]<-savedsims[[j]]$A[[i]][3,17]
      responseSP[i]<-savedsims[[j]]$A[[i]][15, 17]
      responseSF[i]<-savedsims[[j]]$A[[i]][16, 17]
      responseDP[i]<-savedsims[[j]]$A[[i]][5, 17]
      responseLice[i]<-savedsims[[j]]$A[[i]][10, 17]
      responsePRec[i]<- savedsims[[j]]$A[[i]][13, 17]
      modver[i]<-savedsims[[j]]$modver
      TemptoSP[i]<-savedsims[[j]]$TtoSP
    }


  rPRec<-append(rPRec, responsePRec)
  rSP<-append(rSP, responseSP)
  rSF<-append(rSF, responseSF)
  rDP<-append(rDP, responseDP)
  rBF<-append(rBF, responseBF)
  rC<-append(rC, responseC)
  rL<-append(rL, responseL)
  rLice<-append(rLice, responseLice) 
  mv<-append(mv, modver)
  TtoSP<-append(TtoSP, TemptoSP)
sTSP<-append(sTSP, strengthTSP)
sTPR<-append(sTPR, strengthTPR)


}

Combine all the data for SST press perturbations.

SSTresps<-data.frame(unlist(mv), unlist(rL))
colnames(SSTresps)<-c("modver", "Landings")
SSTresps$BigFish<-unlist(rBF)
SSTresps$Catch<-unlist(rC)
SSTresps$ShallowPrey<-unlist(rSP)
SSTresps$SmallFish<-unlist(rSF)
SSTresps$DeepPrey<-unlist(rDP)
SSTresps$PreRecruits<-unlist(rPRec)
SSTresps$Lice<-unlist(rLice)
SSTresps$TtoSP<-unlist(TtoSP)
SSTresps$strengthTSP<-unlist(sTSP)
SSTresps$strengthTPR<-unlist(sTPR)

dim(SSTresps)
## [1] 1280000      12
head(SSTresps)
##   modver      Landings    BigFish         Catch ShallowPrey   SmallFish
## 1  NP_5N -5.567743e+02 -0.3373427 -2.142732e+02 -10.4242266 -13.5343524
## 2  NP_5N  8.170039e-01 -1.4139693 -4.501974e-01  -1.5017343  -5.8845863
## 3  NP_5N  7.833693e-03 -0.2069672 -4.254942e-03   1.4642595  -2.6094508
## 4  NP_5N -2.256247e+00 -2.0706574 -2.137162e+01  -2.7189780  -2.1021427
## 5  NP_5N -4.045399e-01 -0.1520946 -4.219938e-01  -0.4969617  -0.5395686
## 6  NP_5N -8.295007e-01 -2.9448002 -2.165417e+00 -16.1121593  -0.9627455
##     DeepPrey PreRecruits          Lice TtoSP strengthTSP strengthTPR
## 1 0.47505031  -11.643271  1.260652e-15   neg -0.86124457  -0.5398636
## 2 2.02972345   -4.692020  3.508982e-16   neg -0.67097418  -0.8748754
## 3 0.07570563   -4.735434  0.000000e+00   neg -0.02940756  -0.8856475
## 4 0.14613705   -4.989892 -2.823497e-17   neg -0.25066168  -0.6775993
## 5 0.23512084   -2.247589 -1.416418e-15   neg -0.45667532  -0.7669166
## 6 4.24705484    8.161770  0.000000e+00   neg -0.94149840  -0.7914249
#SSTresps2<- SSTresps

We just want the qualitative responses for the purpose of investigating the interplay between edge strength and node responses. The raw response values simulated by QPress in the list A are not meaningful. So we convert them to qualitative responses 1, 0 or -1. (note that QPress provides functions for dealing with this when looking at mean responses across simulations, but here we want to use the data for the strength-response data within each individual simulation. So we need to do it manually.)

Convert responses larger than 0.001 to 1, responses smaller than -0.001 to -1, and everything between to 0.

#Create a new column of qualitative responses for each node.

##Landings
SSTresps$Landings1<-NA
SSTresps$Landings1[SSTresps$Landings<= -1/1000] <- -1
SSTresps$Landings1[SSTresps$Landings> 1/1000] <- 1
SSTresps$Landings1[SSTresps$Landings > -1/1000 & SSTresps$Landings < 1/1000] <- 0
summary(SSTresps$Landings1)       
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -1.0000 -1.0000 -1.0000 -0.2873  1.0000  1.0000
##BigFish
SSTresps$BigFish1<-NA
SSTresps$BigFish1[SSTresps$BigFish<= -1/1000] <- -1
SSTresps$BigFish1[SSTresps$BigFish> 1/1000] <- 1
SSTresps$BigFish1[SSTresps$BigFish > -1/1000 & SSTresps$BigFish < 1/1000] <- 0
summary(SSTresps$BigFish1)  
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -1.0000 -1.0000 -1.0000 -0.4008  1.0000  1.0000
##Catch
SSTresps$Catch1<-NA
SSTresps$Catch1[SSTresps$Catch<= -1/1000] <- -1
SSTresps$Catch1[SSTresps$Catch> 1/1000] <- 1
SSTresps$Catch1[SSTresps$Catch > -1/1000 & SSTresps$Catch < 1/1000] <- 0
summary(SSTresps$Catch1)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -1.0000 -1.0000 -1.0000 -0.4022  1.0000  1.0000
##ShallowPrey
SSTresps$ShallowPrey1<-NA
SSTresps$ShallowPrey1[SSTresps$ShallowPrey<= -1/1000] <- -1
SSTresps$ShallowPrey1[SSTresps$ShallowPrey> 1/1000] <- 1
SSTresps$ShallowPrey1[SSTresps$ShallowPrey > -1/1000 & SSTresps$ShallowPrey < 1/1000] <- 0
summary(SSTresps$ShallowPrey1)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -1.0000 -1.0000  1.0000  0.2746  1.0000  1.0000
##SmallFish
SSTresps$SmallFish1<-NA
SSTresps$SmallFish1[SSTresps$SmallFish<= -1/1000] <- -1
SSTresps$SmallFish1[SSTresps$SmallFish> 1/1000] <- 1
SSTresps$SmallFish1[SSTresps$SmallFish > -1/1000 & SSTresps$SmallFish < 1/1000] <- 0
summary(SSTresps$SmallFish1)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -1.0000 -1.0000 -1.0000 -0.5723 -1.0000  1.0000
##DeepPrey
SSTresps$DeepPrey1<-NA
SSTresps$DeepPrey1[SSTresps$DeepPrey<= -1/1000] <- -1
SSTresps$DeepPrey1[SSTresps$DeepPrey> 1/1000] <- 1
SSTresps$DeepPrey1[SSTresps$DeepPrey > -1/1000 & SSTresps$DeepPrey < 1/1000] <- 0
summary(SSTresps$DeepPrey1)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  -1.000  -1.000   1.000   0.405   1.000   1.000
##PreRecruits
SSTresps$PreRecruits1<-NA
SSTresps$PreRecruits1[SSTresps$PreRecruits<= -1/1000] <- -1
SSTresps$PreRecruits1[SSTresps$PreRecruits> 1/1000] <- 1
SSTresps$PreRecruits1[SSTresps$PreRecruits > -1/1000 & SSTresps$PreRecruits < 1/1000] <- 0
summary(SSTresps$PreRecruits1)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -1.0000 -1.0000 -1.0000 -0.8446 -1.0000  1.0000
##Lice
SSTresps$Lice1<-NA
SSTresps$Lice1[SSTresps$Lice<= -1/1000] <- -1
SSTresps$Lice1[SSTresps$Lice> 1/1000] <- 1
SSTresps$Lice1[SSTresps$Lice > -1/1000 & SSTresps$Lice < 1/1000] <- 0
summary(SSTresps$Lice1)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##       0       0       0       0       0       0

Save the data (qualitative responses to SST press perturbations)

#save(SSTresps, file='Data/SSTresps.Rdata')

This SSTresps data can be used to recreate Fig5a using the script in Fig5a.html