
predict1RFU = function(prot,ano){
	if(!prot%in%names(blood_models)){
		stop(paste0(prot," protein is not part of the Somalogic platform."))
	}
	if(sum(c("GA","Smoker","BMI","Age","Nulliparous")%in%colnames(ano))!=5){
		stop("Annotation table is not in the proper format.")
	}
	colnames(ano)[colnames(ano)=="Smoker"] = "smoker"
	colnames(ano)[colnames(ano)=="Age"] = "age"

	ano$Main_Index=0
	gsplines <- predict(z,ano$GA)
	colnames(gsplines) <- paste("t",colnames(gsplines),sep="")
	ano <- cbind(ano,gsplines)
	predict(blood_models[[prot]],newdata = ano,allow.new.levels=T,re.form=~0)
}



#' predictRFU 
#'
#' Returns the protein abundance expressed as log2 relative fluorescence units (RFU), for the proteins and samples provided by the user.
#'
#' @param protvect 		Vector of proteins for which the log2 RFU are computed.
#' @param ano 			Data frame that contains the sample annotations to be used as predictors by the models, such as gestational age, maternal age, nulliparity, body mass index, and smoking status.
#'
#' @export
predictRFU = function(protvect,ano){
	if(!all(protvect%in%names(blood_models))){
		stop(paste0(paste(setdiff(protvect,names(blood_models)),collapse=", ")," proteins from protvect are not part of the Somalogic platform."))
	}
	rfus = sapply(protvect, function(prot) predict1RFU(prot,ano))
	if(ncol(rfus)==1) {rfus = rfus[,1]}
	rfus
}


#' plotProtTrend 
#'
#' @param prot 				Protein for which the log2 RFU trend is plotted.
#' @param ano 				List of values for the maternal characteristics (Smoker,BMI,Age,Nulliparous) used by the models. The maternal characteristic that has more than one value will be plotted. By default, the Nulliparous trends will be drawn.
#' @param GAlims 			Minimum and maximum for the sequence of gestational age values used as predictor by the models.
#' @param legendgroups 		Character vector for the group names to be used when the group names in the ano list are not sufficiently descriptive.
#' @param legendposition 	Position of the legend in the plot area. By default it is "bottomright".
#' @param palette 			Character vector for the color palette used to plot the trend lines. The default palette has four colors.
#' @param xlim 				x axis minimum and maximum limits
#' @param ylim 				y axis minimum and maximum limits
#'
#' @export
plotProtTrend = function(prot,ano=list(Smoker="0",BMI=26,Age=24,Nulliparous=c("0","1")),
	GAlims=c(8,40),legendgroups=NA,legendposition="bottomright",
	palette=c("blue","red","green","black"),xlim=NA,ylim=NA){

	predictor = names(ano)[sapply(ano,length)!=1]
	if(is.na(legendgroups)){legendgroups=paste(predictor,names(table(ano[predictor])))}
	col	= na.omit(palette[1:length(legendgroups)])
	if(length(ano)!=4|sum(names(ano)%in%c("Smoker","BMI","Age","Nulliparous"))!=4){
		stop("ano argument is not in the proper format.")
	} else if (sum(sapply(ano,length)==1)!=3){
		stop("All but one ano list elements should have a length of 1.")
	} else if (length(ano[[names(ano)[sapply(ano,length)!=1]]])!=length(legendgroups)){
		stop("Predictor used for comparison should have the same number of levels and colors in the legend.")
	}

	ano2 = expand.grid(
		GA = seq(GAlims[1],GAlims[2],0.1),
		smoker = ano$Smoker,
		BMI = ano$BMI,
		age = ano$Age,
		Nulliparous = ano$Nulliparous,
		Main_Index = 0
	)

	if(predictor=="Smoker"){
		predictor="smoker"
	} else if (predictor=="Age"){
		predictor="age"
	}
	levels = names(table(ano2[predictor]))
	for(indl in 1:length(levels)){
		lev = levels[indl]
		nd = ano2[ano2[predictor]==lev,]
		gsplines <- predict(z,nd$GA)
		colnames(gsplines) <- paste("t",colnames(gsplines),sep="")
		nd <- cbind(nd,gsplines)
		expr = predict(blood_models[[prot]],newdata = nd,allow.new.levels=T,re.form=~0)
		if(indl!=1){
			points(nd$GA,expr,col=col[indl],pch=20,cex=0.5)
		} else {
			if(!is.na(xlim)&&!is.na(ylim)) {
				plot(nd$GA,expr,col=col[indl],pch=20,cex=0.5,main=prot,
					xlab="Gestational age (weeks)",ylab="Log2 RFU",xlim=xlim,ylim=ylim)
			} else if(!is.na(xlim)&&is.na(ylim)) {
				plot(nd$GA,expr,col=col[indl],pch=20,cex=0.5,main=prot,
					xlab="Gestational age (weeks)",ylab="Log2 RFU",xlim=xlim)
			} else if(is.na(xlim)&&!is.na(ylim)) {
				plot(nd$GA,expr,col=col[indl],pch=20,cex=0.5,main=prot,
					xlab="Gestational age (weeks)",ylab="Log2 RFU",ylim=ylim)
			} else if(is.na(xlim)&&is.na(ylim)) {
				plot(nd$GA,expr,col=col[indl],pch=20,cex=0.5,main=prot,
					xlab="Gestational age (weeks)",ylab="Log2 RFU")
			}
		}
	}

	if(!is.na(legendposition)){
		legend(legendposition, legend=legendgroups,
			col=col, lty=1, cex=0.8, lwd=2)
	}
}



#' plotRFU 
#' 
#' Plots the protein abundance expressed as log2 relative fluorescence units (RFU) for the protein provided by the user, across all samples.
#'
#' @param prot 			String representing the name of the protein for which the log2 RFU of the samples is plotted.
#' @param protdat 		Data frame that contains the relative fluorescence units (RFU) for the proteins and samples of interest. It has to include the protein prot.
#' @param ano 			Data frame that contains the sample annotations to be used as predictors by the models, such as gestational age, maternal age, nulliparity, body mass index, and smoking status.
#'
#' @export
plotRFU = function(prot,protdat,ano){
	if(!prot%in%rownames(protdat)){
		stop(paste0(prot," protein should be in the rownames of protdat."))
	}
	if(!identical(sort(union(colnames(protdat),rownames(ano))),sort(intersect(colnames(protdat),rownames(ano))))){
		stop("protdat and ano arguments should have the same sample ids.")
	}
	protdat = log2(protdat)
	ano$Y = protdat[prot,rownames(ano)]
	plot(ano$GA,ano$Y,col="grey45",pch=1,lwd=1.6,main=prot,xlab="Gestational age (weeks)",ylab="Log2 RFU")

	galims = range(ano$GA)
	ano = expand.grid(
		GA = seq(galims[1],galims[2],0.1),
		smoker = "0",
		BMI = median(ano$BMI,na.rm=T),
		age = median(ano$Age,na.rm=T),
		Nulliparous = "1",
		Main_Index = 0
	)

	gsplines <- predict(z,ano$GA)
	colnames(gsplines) <- paste("t",colnames(gsplines),sep="")
	ano <- cbind(ano,gsplines)
	expr = predict(blood_models[[prot]],newdata = ano,allow.new.levels=T,re.form=~0)
	points(ano$GA,expr,col="red",lwd=2,type="l")
}



compute1MOM = function(prot,protdat,ano){
	if(sum(c("GA","Smoker","BMI","Age","Nulliparous")%in%colnames(ano))!=5){
		stop("Annotation table is not in the proper format.")
	}
	colnames(ano)[colnames(ano)=="Smoker"] = "smoker"
	colnames(ano)[colnames(ano)=="Age"] = "age"

	protdat = log2(protdat)
	ano2 = ano
	ano2$Y = protdat[prot,rownames(ano2)]
	ano$Main_Index=0
	gsplines <- predict(z,ano$GA)
	colnames(gsplines) <- paste("t",colnames(gsplines),sep="")
	ano <- cbind(ano,gsplines)
	ano2$Y-predict(blood_models[[prot]],newdata = ano,allow.new.levels=T,re.form=~0)
}



#' computeMOM
#' 
#' Returns the log2 multiple of the median (MoM) of protein abundance expressed as relative fluorescence units (RFU), for the proteins and samples provided by the user.
#'
#' @param protvect 		Vector of proteins for which the log2 MOMs of protein abundance are computed.
#' @param protdat 		Data frame that contains the relative fluorescence units (RFU) for the proteins and samples of interest.
#' @param ano 			Data frame that contains the sample annotations to be used as predictors by the models, such as gestational age, maternal age, nulliparity, body mass index, and smoking status.
#'
#' @export
computeMOM = function(protvect=rownames(protdat),protdat,ano){
	if(!all(protvect%in%rownames(protdat))){
		stop("The proteins from protvect should be in the rownames of protdat.")
	}
	if(!identical(sort(union(colnames(protdat),rownames(ano))),sort(intersect(colnames(protdat),rownames(ano))))){
		stop("protdat and ano arguments should have the same sample ids.")
	}
	moms = sapply(protvect, function(prot) compute1MOM(prot,protdat,ano))
	if(ncol(moms)==1) {moms = moms[,1]}
	moms
}



#' plotMOM 
#' 
#' Plots the log2 multiple of the median (MoM) of protein abundance expressed as relative fluorescence units (RFU) for the protein provided by the user, across all samples.
#'
#' @param prot 			String representing the name of the protein for which the log2 MOMs of protein abundance of the samples is plotted.
#' @param protdat 		Data frame that contains the relative fluorescence units (RFU) for the proteins and samples of interest. It has to include the protein prot.
#' @param ano 			Data frame that contains the sample annotations to be used as predictors by the models, such as gestational age, maternal age, nulliparity, body mass index, and smoking status.
#'
#' @export
plotMOM = function(prot,protdat,ano){
	# the conditions on prot ids are enforced by the computeMOM call
	ano$MOM = computeMOM(prot,protdat,ano)
	ano = na.omit(ano)
	plot(ano$GA,ano$MOM,col="grey45",pch=1,lwd=1.6,main=prot,xlab="Gestational age (weeks)",ylab="Log2 MOM")
	abline(h=0,col="black",lty=2,lwd=2.4)
}

