#----------------------------------------------------------
#Détermination des paramètres de l'ajustement

#Décompte des expositions par âge
#Calcul des logits
tauxBruts[is.na(tauxBruts)]=0

resTables[[format(tauxEchantillonnage*100)]]$expo=expo
resTables[[format(tauxEchantillonnage*100)]]$tauxBruts=tauxBruts

logitsBruts=log(tauxBruts/(1-tauxBruts))
logitsBruts[logitsBruts==-Inf]=0

logitsBCAC=log(tauxBCAC/(1-tauxBCAC))
logitsBCAC[logitsBCAC==-Inf]=0

#----------------------------------------------------------
#La régression est effectuée en utilisant comme fonction de perte l'écart absolu entre les qx.

#Fonctions de perte pour l'optimisation
erreur=function(p){
	y=p[1]+p[2]*logitsBCAC[(ancMin+1):ancMax,]
	qAjuste=exp(y)/(1+exp(y))
	sum(abs(expo*(tauxBruts-qAjuste)))
	}

#Optimisation non linéaire
res=constrOptim(c(0,1),erreur,ui=c(0,1),ci=0,control=list(maxit=10^3),method="Nelder-Mead")$par

#----------------------------------------------------------
#On calcule les valeurs ajustées
tauxAjustes=exp(res[1]+res[2]*logitsBCAC[(ancMin+1):ancMax,])/(1+exp(res[1]+res[2]*logitsBCAC[(ancMin+1):ancMax,]))
resTables[[format(tauxEchantillonnage*100)]][[modeleBrass]]$tauxAjustes=tauxAjustes

#----------------------------------------------------------
#Variante W-H
library(RColorBrewer)
source(paste(dossierOutils,"\\WH_2.R",sep=""))

tauxBruts35=tauxBruts[1:35,]
expo35=expo[1:35,]

tauxAjustes35=W_H(tauxBruts35,expo35,2,2,1000,1000)
tauxAjustes=rbind(tauxAjustes35,pmin(1,c(1:40)))
tauxAjustes=pmax(tauxAjustes,0)
resTables[[format(tauxEchantillonnage*100)]][[modeleWH]]$tauxAjustes=tauxAjustes

Ecart=tauxBruts35/tauxAjustes35

q=ncol(tauxBruts35)
p=nrow(tauxBruts35)
x=matrix(1:p,p,1)
y=matrix(1:q,q,1)

par(mfrow=c(1,2))
sauveFigure()
persp(x,y,tauxBruts35,theta=-30,col=brewer.pal(9,"RdYlGn"),xlab="Age",ylab="Année",main="Taux bruts",adj=0.5,font=2)
persp(x,y,tauxAjustes35,theta=-30,col=brewer.pal(9,"RdYlGn"),xlab="Age",ylab="Année",main="Taux lissés",adj=0.5,font=2)
dev.off()
#----------------------------------------------------------
#Comparaison des nombres de sorties modélisés et observés globalement
sortiesObservees=sum(expo*tauxBruts)
sortiesTheoriques=sum(expo*tauxAjustes)


