##################################################################################################################
# Calcula la media y la desviación típica de una serie de columnas de un dataframe, así como sus errores estandar
# sdata: dataframe con los datos
# cmean: vector de índices o nombres de columnas sobre las cuales queremos realizar los cálculos
# wght: índice o nombre de la columna que contiene el peso de la fila
# brr: vector con los índices o nombres de columna con los pesos replicados
# Devuelve una matriz con cuatro filas y una columna para cada columna indicada en cmean
# la primera fila contiene las medias, la segunda el error estandar de la media, la tercera la desviación típica 
# y la cuarta el error estandar de la desviación típica
##################################################################################################################
# Calculates the mean and standard deviation of some columns of a dataframe, along with their standard errors
# sdata: dataframe with the data
# cmean: vector with the indices or names of the columns on which perform the calculations
# wght: index or name of the column with the weight of the rows
# brr: vector with the indices or names of the columns with the replicate weights
# Returns a matrix with four rows and a column for each of the columns in the cmean parameter
# the firs row is for the mean, the second for the mean standard error, the third for the standard deviation
# and the last one for the standard deviation standard error

wght_meansd<-function(sdata,cmean,wght,brr) {
	mmeans<-matrix(ncol=length(cmean),nrow=4);
	mmeans[,]<-0;
	if (is.numeric(cmean)) {
		colnames(mmeans)<-names(sdata)[cmean];
	}
	else {
		colnames(mmeans)<-cmean;
	}	
	rownames(mmeans)<-c("MEAN","SE-MEAN","STDEV","SE-STDEV");
	swght<-sum(sdata[,wght]);
	for (i in 1:length(cmean)) {
		mmeans[1,i]<-sum(sdata[,wght]*sdata[,cmean[i]])/swght;
		mmeans[3,i]<-sqrt((sum(sdata[,wght]*(sdata[,cmean[i]]^2))/swght)-mmeans[1,i]^2);
		for (j in 1:length(brr)) {
			sbrr<-sum(sdata[,brr[j]]);
			mmbrj<-sum(sdata[,brr[j]]*sdata[,cmean[i]])/sbrr;
			mmeans[2,i]<-mmeans[2,i] + (mmbrj - mmeans[1,i])^2;
			mmeans[4,i]<-mmeans[4,i] + (sqrt((sum(sdata[,brr[j]]*(sdata[,cmean[i]]^2))/sbrr)-mmbrj^2) - mmeans[3,i])^2;
		}		
	}
	mmeans[2,]<-sqrt((mmeans[2,] * 4) / length(brr));
	mmeans[4,]<-sqrt((mmeans[4,] * 4) / length(brr));
	return(mmeans);
}

##################################################################################################################
# Calcula la media y la desviación típica por país de una serie de columnas de un dataframe, así como sus errores estandar
# sdata: dataframe con los datos
# cmean: vector de índices o nombres de columnas sobre las cuales queremos realizar los cálculos
# cnt: índice o nombre de la columna con el país
# wght: índice o nombre de la columna que contiene el peso de la fila
# brr: vector con los índices o nombres de columna con los pesos replicados
# Devuelve una lista de matrices, una por país, cada una con con cuatro filas y una columna para cada columna indicada en cmean
# la primera fila contiene las medias, la segunda el error estandar de la media, la tercera la desviación típica 
# y la cuarta el error estandar de la desviación típica
##################################################################################################################
# Calculates the mean and standard deviation by country of some columns of a dataframe, along with their standard errors
# sdata: dataframe with the data
# cmean: vector with the indices or names of the columns on which perform the calculations
# cnt: index or name of the column with the country
# wght: index or name of the column with the weight of the rows
# brr: vector with the indices or names of the columns with the replicate weights
# Returns a list of matrixex, one for country, each one with four rows and a column for each of the columns in the cmean parameter
# the firs row is for the mean, the second for the mean standard error, the third for the standard deviation
# and the last one for the standard deviation standard error

wght_meansdbycnt<-function(sdata,cmean,cnt,wght,brr) {
	lcntrs<-vector('list',length(levels(as.factor(sdata[,cnt]))));
	for (p in 1:length(levels(as.factor(sdata[,cnt])))) {
		names(lcntrs)[p]<-levels(as.factor(sdata[,cnt]))[p];
	}
	for (p in 1:length(levels(as.factor(sdata[,cnt])))) {
		mmeans<-matrix(ncol=length(cmean),nrow=4);
		mmeans[,]<-0;
		if (is.numeric(cmean)) {
			colnames(mmeans)<-names(sdata)[cmean];
		}
		else {
			colnames(mmeans)<-cmean;
		}		
		rownames(mmeans)<-c("MEAN","SE-MEAN","STDEV","SE-STDEV");
		rcnt <- sdata[,cnt]==levels(as.factor(sdata[,cnt]))[p];
		swght<-sum(sdata[rcnt,wght]);
		for (i in 1:length(cmean)) {			
			mmeans[1,i]<-sum(sdata[rcnt,wght]*sdata[rcnt,cmean[i]])/swght;
			mmeans[3,i]<-sqrt((sum(sdata[rcnt,wght]*(sdata[rcnt,cmean[i]]^2))/swght)-mmeans[1,i]^2);
			for (j in 1:length(brr)) {
				sbrr<-sum(sdata[rcnt,brr[j]]);
				mmbrj<-sum(sdata[rcnt,brr[j]]*sdata[rcnt,cmean[i]])/sbrr;
				mmeans[2,i]<-mmeans[2,i] + (mmbrj - mmeans[1,i])^2;
				mmeans[4,i]<-mmeans[4,i] + (sqrt((sum(sdata[rcnt,brr[j]]*(sdata[rcnt,cmean[i]]^2))/sbrr)-mmbrj^2) - mmeans[3,i])^2;
			}		
		}
		mmeans[2,]<-sqrt((mmeans[2,] * 4) / length(brr));
		mmeans[4,]<-sqrt((mmeans[4,] * 4) / length(brr));
		lcntrs[[p]]<-mmeans;
	}
	return(lcntrs);
}

##################################################################################################################
# Calcula la media y la desviación típica de una serie de columnas de un dataframe, así como sus errores estandar
# los cálculos para cada columna se desagregan por una serie de factores
# sdata: dataframe con los datos
# cmean: vector de índices o nombres de columnas sobre las cuales queremos realizar los cálculos
# cfact: vector con los índices o nombres de columnas con los factores para desagregar los cálculos
# wght: índice o nombre de la columna que contiene el peso de la fila
# brr: vector con los índices o nombres de columna con los pesos replicados
# Devuelve una matriz con cuatro filas y una columna para cada combinación columna / valor factor indicada en cmean y cfact
# la primera fila contiene las medias, la segunda el error estandar de la media, la tercera la desviación típica 
# y la cuarta el error estandar de la desviación típica
##################################################################################################################
# Calculates the mean and standard deviation of some columns of a dataframe, along with their standard errors
# the calculations for each column will be disaggregated by some factors
# sdata: dataframe with the data
# cmean: vector with the indices or names of the columns on which perform the calculations
# cfact: vector with the indices or names of the columns with the factors to disaggregate the calculations
# wght: index or name of the column with the weight of the rows
# brr: vector with the indices or names of the columns with the replicate weights
# Returns a matrix with four rows and a column for each of the combination column / factor level in the cmean / cfact parameters
# the firs row is for the mean, the second for the mean standard error, the third for the standard deviation
# and the last one for the standard deviation standard error

wght_meansdfact<-function(sdata,cmean,cfact,wght,brr) {
	nc<-0;
	for (i in 1:length(cfact)) {
		nc <- nc + length(levels(as.factor(sdata[,cfact[i]])));
	}
	mmeans<-matrix(ncol=nc,nrow=4 * length(cmean));
	mmeans[,]<-0;
	cn<-c();
	for (i in 1:length(cfact)) {
		for (j in 1:length(levels(as.factor(sdata[,cfact[i]])))) {
			cn<-c(cn, paste(names(sdata)[cfact[i]],levels(as.factor(sdata[,cfact[i]]))[j],sep="-"));
		}
	}
	colnames(mmeans)<-cn;
	rn<-c();
	for (i in 1:length(cmean)) {
		rn<-c(rn, paste(names(sdata)[cmean[i]],"MEAN",sep="-"));
		rn<-c(rn, paste(names(sdata)[cmean[i]],"MEAN","SE",sep="-"));
		rn<-c(rn, paste(names(sdata)[cmean[i]],"STDEV",sep="-"));
		rn<-c(rn, paste(names(sdata)[cmean[i]],"STDEV","SE",sep="-"));
	}	
	rownames(mmeans)<-rn;	
	ic<-1;
	for(f in 1:length(cfact)) {
		for (l in 1:length(levels(as.factor(sdata[,cfact[f]])))) {
			ir<-1;
			rfact<-sdata[,cfact[f]]==levels(as.factor(sdata[,cfact[f]]))[l];
			swght<-sum(sdata[rfact,wght]);
			for (i in 1:length(cmean)) {				
				mmeans[ir,ic]<-sum(sdata[rfact,wght]*sdata[rfact,cmean[i]])/swght;
				mmeans[ir+2,ic]<-sqrt((sum(sdata[rfact,wght]*(sdata[rfact,cmean[i]]^2))/swght)-mmeans[ir,ic]^2);
				for (j in 1:length(brr)) {
					sbrr<-sum(sdata[rfact,brr[j]]);
					mmbrj<-sum(sdata[rfact,brr[j]]*sdata[rfact,cmean[i]])/sbrr;
					mmeans[ir+1,ic]<-mmeans[ir+1,ic] + (mmbrj - mmeans[ir,ic])^2;
					mmeans[ir+3,ic]<-mmeans[ir+3,ic] + (sqrt((sum(sdata[rfact,brr[j]]*(sdata[rfact,cmean[i]]^2))/sbrr)-mmbrj^2) - mmeans[ir+2,ic])^2;
				}		
				mmeans[ir+1,ic]<-sqrt((mmeans[ir+1,ic] * 4) / length(brr));
				mmeans[ir+3,ic]<-sqrt((mmeans[ir+3,ic] * 4) / length(brr));
				ir<-ir+4;
			}
			ic<-ic + 1;
		}
	}
	return(mmeans);
}