misc/Prod_indicateur_tc.R

In cesco-lab/Vigie-Chiro_scripts: Predict Distribution Of Chiropter

library(data.table)
#A EDITER en fonction de la localisation des fichiers .tc
#TC_path="C:/Users/Yves Bas/Documents/txt"
#get arguments from the command line
args <- commandArgs(trailingOnly = TRUE)
#uncomment the following line if you prefer to do not use R in command line
if(length(args)==0){
  args="D:/DataTestVeolia/nuit_du_23-06-2019_au_24-06-2019"
  #q()
}
TC_path=args[1]
Ref_coef=fread("RefPF_Coef (5).csv")
DataCalibr=fread("DataCalibr.csv")

alogit <- function(x) 
{
  exp(x)/(1+exp(x))
}

#get the .tc files list
idlist=list.files(TC_path,pattern=".tc$",full.names=T,recursive=T)
print(length(idlist))
my.data <- list()
for(f in 1:length(idlist)) {
  my.data[[f]] <- fread(idlist[[f]])
}
Idtot=as.data.frame(rbindlist(my.data))

SpMax<-max.col(Idtot[,2:(ncol(Idtot)-6)],ties.method = "first")
ScoreMax=apply(Idtot[,2:(ncol(Idtot)-6)],MARGIN=1,FUN=max)
SpMax2=as.data.frame(cbind(Filename=as.character(Idtot[,1]),Id=colnames(Idtot)[SpMax+1],numsp=SpMax,IndConf=ScoreMax))


SpMax3=merge(SpMax2,Ref_coef,by.x="Id",by.y="Code")
SpMax3$IndConf=as.numeric(as.character(SpMax3$IndConf))

test=match(Ref_coef$Code,levels(as.factor(SpMax3$Id)))

Sub_Ref_coef=subset(Ref_coef,is.na(test)==F)
if(nrow(SpMax3)==0)
{
  print("ATTENTION : aucunes chauves-souris contactées !")
  Indicateurs=c(0,0,0)
}else{

  test2=match(levels(as.factor(SpMax3$Id)),Sub_Ref_coef$Code)
  
fwrite(Sub_Ref_coef,"/log0.csv",sep=";")
  
if (exists("DataInd")){rm(DataInd)}
ListSp=vector()  
AbBrute=vector()
for (i in 1:nlevels(as.factor(SpMax3$Id)))
{
  Datasub=subset(SpMax3,(SpMax3$Id==levels(as.factor(SpMax3$Id))[i])
                 &(as.numeric(as.character(SpMax3$IndConf))>SpMax3$Seuil))
  Ind_S=log10(nrow(Datasub)+1)/Sub_Ref_coef$Moy[test2[i]]*Sub_Ref_coef$SSIZ[test2[i]]
  Ind_L=-log10(nrow(Datasub)+1)/Sub_Ref_coef$Moy[test2[i]]*Sub_Ref_coef$ReponseLum[test2[i]]
  Ind_A=log10(nrow(Datasub)+1)/Sub_Ref_coef$Moy[test2[i]]*Sub_Ref_coef$ReponseArb[test2[i]]
  Ind_E=log10(nrow(Datasub)+1)/Sub_Ref_coef$Moy[test2[i]]*Sub_Ref_coef$ReponseEau[test2[i]]
  AbStand=log10(nrow(Datasub)+1)/Sub_Ref_coef$Moy[test2[i]]
  ListSp=c(ListSp,levels(as.factor(SpMax3$Id))[i])
  AbBrute=c(AbBrute,nrow(Datasub))
  if (exists("DataInd")){DataInd=rbind(DataInd,cbind(Ind_S,Ind_L,Ind_A,Ind_E,AbStand))}else{DataInd=cbind(Ind_S,Ind_L,Ind_A,Ind_E,AbStand)}
  print(paste(levels(as.factor(SpMax3$Id))[i],nrow(Datasub),Ind_S,Ind_L,Ind_A,Ind_E,AbStand))
}

fwrite(DataInd,"/log1.csv",sep=";")


DataInd[is.na(DataInd)]=0
DataInd=as.data.frame(DataInd)
DataInd=subset(DataInd,(DataInd$AbStand!=Inf)&(DataInd$AbStand>0))

fwrite(DataInd,"/log2.csv",sep=";")
fwrite(DataCalibr,"/log3.csv",sep=";")


  Indicateur_Spec=alogit(((sum(DataInd[,1])/sum(DataInd[,5]))-DataCalibr[1,1])/DataCalibr[1,2])*20
  Indicateur_Lum=alogit((log(sum(DataInd[,2])+0.1)-DataCalibr[2,1])/DataCalibr[2,2])*20
  Indicateur_Arb=alogit((log(sum(DataInd[,3])+0.5)-DataCalibr[3,1])/DataCalibr[3,2])*20
  Indicateur_Eau=alogit((log(sum(DataInd[,3])+2)-DataCalibr[4,1])/DataCalibr[4,2])*20
  
  
BilanEspeces0=as.data.frame(cbind(ListSp,AbBrute)) 
BilanEspeces1=merge(BilanEspeces0,Ref_coef,by.x="ListSp",by.y="Code",all.y=T)
BilanEspeces1$AbBrute=as.character(BilanEspeces1$AbBrute)
BilanEspeces1$AbBrute[is.na(BilanEspeces1$AbBrute)]=0
BilanEspeces=as.data.frame(cbind(Groupe=BilanEspeces1$Groupe,Espece=BilanEspeces1$Espece,Nb_Contacts=BilanEspeces1$AbBrute))
BilanEspeces=BilanEspeces[order(as.character(BilanEspeces$Groupe),BilanEspeces$Espece),]

Indicateurs=c(Indicateur_Spec,Indicateur_Lum,Indicateur_Arb,Indicateur_Eau)

#A EDITER si l'on veut changer la destination des données
write.csv(Indicateurs,paste(TC_path, "/Indicateurs.csv", "", sep=""),row.names=F)
write.csv(BilanEspeces,paste(TC_path, "/BilanEspeces.csv", "", sep=""),row.names=F)
}