Messagepar hénin virginie » 30 Mai 2006, 07:57
ben j'ai pas vraiment l'habitude du C mais i tu pouvais jeter un coup d'oeil à mon programme et me dire si tu vois des amméliorations merci
beaucoup
voilà mon code alors ::
new_ponderation=function(chemin1,chemin2,chemin3)
{
#chemin1 : chemin permettant d'accéder au premier tableau de sortie de la BAO
#chemin2 : chemin permettant d'accéder aux tableaux rbao* de sortie de la BAO
#contient les données restructurées à l'aide de perl
#début de la fonction de vraisemblance
cat("Sur combien de campagnes porte l'analyse\n")
nbcamp<-scan(quiet=T,"",numeric(),1)
Nb_var=length(xvar)
cat("Sur combien d'échantillons calcule t'on la vraisemblance \n")
Nb_ech<-scan(quiet=T,"",numeric(),1)
tab_b=read.table(chemin1,sep=" ")
tab_b=cbind(as.factor(tab_b[,1]),tab_b[4:5],tab_b[7])
names(tab_b)=c("N_camp","N_var","nom_var","Kij")
Nb_param=length(xname)
print ("debut programme")
print(date())
sigma=matrix(0,Nb_var,Nb_ech)
somme_obs_predi=matrix(0,Nb_var,Nb_ech)
poids_simu=matrix(0,1,Nb_ech)
poids_log=matrix(0,1,Nb_ech)
#calcul du nombre d'observations par type de mesure
nt=0
for (i in 0:Nb_var-1)
{
nt[i+1]= sum(tab_b[tab_b$N_var==(i),4])
}
tab_nt=cbind(as.vector(tab_b$nom_var[1:Nb_var]),nt)
tab_var=list()
for (i in 1:Nb_var)
{
tab_var[[i]]=as.data.frame(matrix(data =scan(file=paste(chemin3,xvar[i],".txt",sep=""),what = numeric(0)),nrow=as.numeric(tab_nt[i,2])*Nb_ech, ncol = 3, byrow = TRUE))
}
for (e in 1:Nb_ech)
{
#calcul de la pondération pour un échantillon (correspond à une ligne du fichier .rbao
#calcul des sommes observés - prédits au carré
somme_obs_pred=matrix(0,1,Nb_var)
#tirage des 1/sigma^2
sigmaj=matrix(0,1,Nb_var)
poids_simu[e]=1
poids_log[e]=0
for (i in 1:Nb_var)
{ if (dim(tab_var[[i]])[[1]]!=0)
{
tab_var_inter=tab_var[[i]]
somme_obs_pred[i]=sum((tab_var_inter[tab_var_inter$V1 == e,][,2]- tab_var_inter[tab_var_inter$V1 == e,][,3])^2)
sigmaj[i]=rgamma(1,shape=as.numeric(tab_nt[i,2]),scale=2/as.numeric(somme_obs_pred[i]))
poids_log[e]=poids_log[e]+((as.numeric(tab_nt[i,2])/2+1)*log(2/as.numeric(somme_obs_pred[i])))
poids_simu[e]=poids_simu[e]*((2/as.numeric(somme_obs_pred[i]))^(as.numeric(tab_nt[i,2])/2+1))
}
else
{
somme_obs_pred[i]="NA"
sigmaj[i]="NA"
}
}
sigma[,e]=sigmaj
somme_obs_predi[,e]=somme_obs_pred
print (e)
}
param=as.data.frame(matrix(data =scan(file=paste(chemin3,"para.txt",sep=""),what = numeric(0)),nrow=Nb_ech, ncol = Nb_param, byrow = TRUE))
resul=list(poids_simu,param,sigma,somme_obs_predi,poids_log)
names(resul)=c("ponderation","param","sigma","somm_erreur","ponderation_log")
print ("fin programme")
print(date())
return(resul)
}
Cordialement