Commit 289e1716 authored by Henry Gerard's avatar Henry Gerard
Browse files

avec gcc4

parents
#!/bin/bash
#SBATCH --job-name=exoBSC1
#SBATCH -p q-128Go
#SBATCH -o /scratch/frederic.gosselin/%j.out
#SBATCH -e /scratch/frederic.gosselin/%j.err
#SBATCH --mail-type=ALL
#SBATCH --mail-user=frederic.gosselin@irstea.fr
module load gcc/5.2.0
cd /nfs/home/frederic.gosselin/save/BSC
/cm/shared/apps/R/R-3.2.4/bin/R CMD BATCH --no-save exoBSC1.txt
#mme chose qu'avant sauf que control est comme DE616
reference.letter<-"E"
#nom d'utilisateur sr station de calcul Clermont: peut tre "zilliox" ou "gosselin"
utilisateur<-"frederic.gosselin"
#NOTE: il faut pour que cela fonctionne sur la sttaion de calcul:
# sous le "repertoire_commun": la dernire version de dist_packV01.c (mais a priori ne devrait pas changer)
#sous le rpertoire "IFN": la dernire version de "DREAM_modFG.RData", la dernire version de "jeu_entier.csv", et le fichier de commandes "exoIFN.txt"
#if (Sys.info()["sysname"]=="Windows") library(lme4)
library(gdata)
library(coda)
#library(rjags)
#library(R2WinBUGS)
library(lme4)
library(AICcmodavg)
#library(Hmisc)
#if (floor(as.double(version$minor))<=13&floor(as.double(version$major))<=2) {library(Design)} else { library(rms) }
library(rstan)
library(nlme)
library(VGAM)
rstan_options(auto_write = TRUE)
options(mc.cores = parallel::detectCores())
load("BSC_Data_functions.RData")
#################### lanc":
##"
runStan<-T
runDIC<-F
DIC.extinct.fourthsample.supp.muchlonger<-matrix(NA,ncol=5,nrow=8)
dimnames(DIC.extinct.fourthsample.supp.muchlonger)<-list(c(
"POPd","POPdl","POPdsq","GDP1900a","GDP1900al",
"SEAL.2009sq","SEAL.2009l","SEAL.2009"),c("DIC","DICvar","DICmin","IC","PDvar"))
DICmarginal.extinct.fourthsample.supp.muchlonger<-matrix(NA,ncol=5,nrow=8)
dimnames(DICmarginal.extinct.fourthsample.supp.muchlonger)<-list(c(
"POPd","POPdl","POPdsq","GDP1900a","GDP1900al",
"SEAL.2009sq","SEAL.2009l","SEAL.2009"),c("DIC","DICvar","DICmin","IC","PDvar"))
DICmarginal1000r.extinct.fourthsample.supp.muchlonger<-matrix(NA,ncol=5,nrow=8)
dimnames(DICmarginal1000r.extinct.fourthsample.supp.muchlonger)<-list(c(
"POPd","POPdl","POPdsq","GDP1900a","GDP1900al",
"SEAL.2009sq","SEAL.2009l","SEAL.2009"),c("DIC","DICvar","DICmin","IC","PDvar"))
model.delta<-glm(cbind(Extinct,Native-Extinct)~1,family=binomial,data=data.unlisted.sp,
subset=!is.element(COU,c("CYP","LUX","MAL"))&!is.na(SEAL.2009+INC.SEAL.CA+PIB1900aD)&!is.na(Native)&!is.na(Extinct))
model.eta<-glm(cbind(Extinct,Native-Extinct)~Taxon+I(Endemic/Native/0.05):Taxon-1,family=binomial,data=data.unlisted.sp,
subset=!is.element(COU,c("CYP","LUX","MAL"))&!is.na(SEAL.2009+INC.SEAL.CA+PIB1900aD)&!is.na(Native)&!is.na(Extinct))
model.bidon<-glm(Extinct~Native-1,family=gaussian,data=data.unlisted.sp,
subset=!is.element(COU,c("CYP","LUX","MAL"))&!is.na(SEAL.2009+INC.SEAL.CA+PIB1900aD)&!is.na(Native)&!is.na(Extinct))
model.bidonbis<-glm(Native~Extinct-1,family=gaussian,data=data.unlisted.sp,
subset=!is.element(COU,c("CYP","LUX","MAL"))&!is.na(SEAL.2009+INC.SEAL.CA+PIB1900aD)&!is.na(Native)&!is.na(Extinct))
model.COU<-glm(cbind(Extinct,Native-Extinct)~COU-1,family=binomial,data=data.unlisted.sp,
subset=!is.element(COU,c("CYP","LUX","MAL"))&!is.na(SEAL.2009+INC.SEAL.CA+PIB1900aD)&!is.na(Native)&!is.na(Extinct))
model.Ptype<-glm(cbind(Extinct,Native-Extinct)~1,family=binomial,data=data.unlisted.sp,
subset=!is.element(COU,c("CYP","LUX","MAL"))&!is.na(SEAL.2009+INC.SEAL.CA+PIB1900aD)&!is.na(Native)&!is.na(Extinct))
model.sdinner<-glm(cbind(Extinct,Native-Extinct)~Taxon-1,family=binomial,data=data.unlisted.sp,
subset=!is.element(COU,c("CYP","LUX","MAL"))&!is.na(SEAL.2009+INC.SEAL.CA+PIB1900aD)&!is.na(Native)&!is.na(Extinct))
win.data <-list(Y = as.vector(model.matrix(model.bidonbis)),
N = as.vector(model.matrix(model.bidon)),
matrixeta=model.matrix(model.eta),
matrixdelta=model.matrix(model.delta),
matrixsdinner=model.matrix(model.sdinner),
Ncases=length( (model.delta)$y),
Ndelta=dim(model.matrix(model.delta))[2],
Neta=dim(model.matrix(model.eta))[2],
Nsdinner=dim(model.matrix(model.sdinner))[2],
COU=apply(model.matrix(model.COU),1,function(x){which(x>0)}),
Ptype=apply(model.matrix(model.COU),1,function(x){which(x>0)}),
NCOU=length(unique(apply(model.matrix(model.COU),1,function(x){which(x>0)}))))
nc <- 3 ;
ni <- 350000;
nb <- 20000 ;
nt <- 30;
params <- c("pareta","parsdinner","sdCOU","easCOU","log_lik")
if (runStan)
{
set.seed(29)
modExtinct_betabinomheteroscpriorsdlowrest_Null_endemicTaxon_fourthsample_bis <- stan(model_code = stanmodel.heteroscbetabinom.priorsdlowrestricted, model_name = "Stan_binom", data = win.data,init = inits.betabinomquad.extinct,pars = params,iter = ni, chains = nc, warmup=nb, thin = nt, save_dso = TRUE,seed = 28,verbose = TRUE)
#
save(modExtinct_betabinomheteroscpriorsdlowrest_Null_endemicTaxon_fourthsample_bis,file="modExtinct_betabinomheteroscpriorsdlowrest_Null_endemicTaxon_fourthsample_bis.RData")
rm(modExtinct_betabinomheteroscpriorsdlowrest_Null_endemicTaxon_fourthsample_bis)
}
if (runDIC)
{attach("modExtinct_betabinomheteroscpriorsdlowrest_Null_endemicTaxon_fourthsample_bis.RData",pos=2)
DIC.extinct.fourthsample.supp.muchlonger<-rbind(DIC.extinct.fourthsample.supp.muchlonger,Null_endemicTaxon=DIC.heteroscbetabinom(modExtinct_betabinomheteroscpriorsdlowrest_Null_endemicTaxon_fourthsample_bis,win.data=win.data))
DICmarginal.extinct.fourthsample.supp.muchlonger<-rbind(DICmarginal.extinct.fourthsample.supp.muchlonger,Null_endemicTaxon=DICmarginal.stanmodel.heteroscbetabinom(modExtinct_betabinomheteroscpriorsdlowrest_Null_endemicTaxon_fourthsample_bis,win.data=win.data))
DICmarginal1000r.extinct.fourthsample.supp.muchlonger<-rbind(DICmarginal1000r.extinct.fourthsample.supp.muchlonger,Null_endemicTaxon=DICmarginal.stanmodel.heteroscbetabinom(modExtinct_betabinomheteroscpriorsdlowrest_Null_endemicTaxon_fourthsample_bis,win.data=win.data,rep.inner=1000))
detach(2)
#load("E:\\Dossier Frederic\\Articles\\BSC_Callois\\DATA_ANALYSIS\\BSC_Extinct_fourthsample_supp_muchlonger.RData")
save.image("E:\\Dossier Frederic\\Articles\\BSC_Callois\\DATA_ANALYSIS\\BSC_Extinct_fourthsample_supp_muchlonger.RData")
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment