#These are the comands to execute in R to make a SRC sensitivity analysis of the different models
#Sets of parameters are alreay defined for each model. Be sure to have the files of parameters definitions in the correct directory. Format : XXX_Parameters.csv.
#Type of files :
#X            min max base unit
#parameter 1   x   x   x    u
#...
#parameter p   x   x   x    u

#Package "sensitivity" needed (available on CRAN)
library(sensitivity,lib.loc="/home/paulb/R/x86_64-redhat-linux-gnu-library/2.15");

#Parameters of the program you can change
########################################################################
XXX="WCT" #WCT, MAC, AM or ABP. Select the model to be tested
Dist="Norm99" #Norm95 or Norm99 or Unif. Select the distribution you want for the parameters. Norm99 is the normal distribution with less than 1% chance of being out of the defined interval.
Directory_Entrie<-"/data/paulb/Data/Parameters/Entries_Sensitivity/" #Directory where the parameters are defined
Directory_Output<-"/data/paulb/Results/Communication/" #Directory where you want your output to be stocked. Has to match with the entry  directory of Scilab.
Directory_Results<-"/data/paulb/Results/Results/" #Directory where you want your results to be stocked.

n <- 1000 # total cost of (p + 2) * n model evaluations
nb <- 10000 #nboot
########################################################################

#Reading the Parameters
Pa<-read.csv(paste(Directory_Entrie,XXX,"_Parameters.csv",sep=""),sep=";"); 
param_min<-Pa$min # Definition of the parameters
param_max<-Pa$max
p<-nrow(Pa);

source("/data/paulb/R/Functions/Random_Param.R") #see comments in file "Random_Param.r"
X<-data.frame(Random_Param(n,param_base))
write.table(X,file=paste(Directory_Output,XXX,"SRC_ParamToBeTested_",n,"_",nb,".csv",sep=""),quote=FALSE,sep=",",row.name=FALSE,col.name=FALSE)
save.image(file=paste("//data/paulb/R/",XXX,"_SRC_ParamToBeTested_",n,"_",nb,"_UNFINISHED.RData",sep=""))

###############################################################
#Execute Scilab Functions infected_XXX.sci and Compute_Sol.sce
###############################################################

#finish the analisis by taking the outputs generated by Scilab
t=9000-2
Y <- as.matrix(read.csv("/Users/paulb/R/LinearRegression/ModelResults.csv",header = FALSE,sep = ","))
All_SRC_squared<-NULL;
All_SRC<-NULL;
CI_squared<-NULL;
CI_squared_col<-NULL;
a<-NULL; b<-NULL;
#WCT:
names(X)<-c("muh", "h","l","alpha","k","v","nuh","r","f_u","g_u","f_n","g_n","x")
for (k in 1:t) {
	x <- src(X, Y[k,], nboot = 100)#compute the SRC indices
	All_SRC<-rbind(All_SRC,x$SRC);
}
write.table(CI_squared,file="WCT_All_SRC_CI_Squared.csv",quote=FALSE,sep=",")#Save the SRC coeficients
write.table(All_SRC,file="WCT_All_SRC.csv",quote=FALSE,sep=",")#Save the SRC coeficients
Rsquare=NULL;
for (i in 1:t){
	Rsquare<-c(Rsquare, sum(All_SRC_squared[,i])); #Add the sum of the SRC squared (equal to R square)
	}
All_SRC_squared<-rbind(All_SRC_squared,Rsquare);
All_SRC_squared<-as.data.frame(All_SRC_squared);
rownames(All_SRC_squared)<-c("muh", "h","l","alpha","k","v","nuh","r","f_u","g_u","f_n","g_n","x","R_square");
write.table(All_SRC_squared,file="WCT_All_SRC_squared_bis.csv",quote=FALSE,sep=",")#Save the squared SRC coeficients

#Plot the results
source(file="/Users/paulb/R/Functions/PlotAllSquared.r")
names<-c("muh", "h","l","alpha","k","v","nuh","r","f_u","g_u","f_n","g_n","x","R_square") #Choose here the parameters to be plotted
xlim<-c(0,t) #Choose the range of times for your plot
ylim<-c(0,1) #Choose the range for y
plotAllSquared(All_SRC,q,names,xlim,ylim)



#MAC:
#names(X)<-c("alpha","b","HD","WN","x","f_U","g_U","l","nuh","muh","n")
#lin<-lm(Y ~ alpha+b+HD+WN+x+f_U+g_U+l+nuh+muh+n,X);