/*  
 *  src/org/metaqtl/MetaQtlResult.java
 * 
 *  $Author: veyrieras $, $Date: 2006-10-17 15:41:54 $, $Version$
 *
 *  MetaQTL : a Java package for QTL Meta-Analysis
 * 
 *  Copyright (C) 2003-2006  Jean-Baptiste Veyrieras, INRA, France.
 *
 *  This program is free software; you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation; either version 2 of the License, or
 *  (at your option) any later version.
 *
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU General Public License for more details.
 *
 *  You should have received a copy of the GNU General Public License
 *  along with this program; if not, write to the Free Software
 *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA * 
 */
package org.metaqtl;

import java.util.Iterator;
import java.util.List;

import org.metaqtl.adapter.QtlAdapter;

/**
 * A <code>MetaQtlResult</code> is used to store the 
 * results of a QTL meta analysis which has been performed
 * for a given trait group on a given chromosome.
 */
public class MetaQtlResult {
	/**
	 * The name of the trait.
	 */
	public String trait = null;
	/**
	 * The number of QTL
	 */
	public int nqtl = 0;
	/**
	 * The number of meta qtl clusterings.  
	 */
	public int nmqc = 0;
	/**
	 * The indices of QTL which are implied
	 * in the analysis.
	 */
	public int[] qtlIdx = null;
	/**
	 * The QTL positions.
	 */
	public double[] x = null;
	/**
	 * The QTL position standard deviations.
	 */
	public double[] sd = null;
	/**
	 * The meta qtl clusterings.
	 */
	public EMResult[] clusterings = null;
	/**
	 * The qtl tree
	 */
	public Tree tree = null;
	
	
	public MetaQtlResult() {};
	/**
	 * @param qtlNumber
	 */
	public MetaQtlResult(String trait, int qtlNumber) {
		nqtl   = qtlNumber;
		qtlIdx = new int[nqtl];
		x      = new double[nqtl];
		sd     = new double[nqtl];
		this.trait = trait;
	}
	/**
	 * Returns the Qtl clustering result with <code>k</code>
	 * clusters. 
	 * @param k
	 * @return
	 */
	public EMResult getClustering(int k) {
		if (clusterings==null || clusterings.length < k+1 || k<0) return null;
		if (k<0) return clusterings[0];
		if (k>nmqc) return clusterings[nmqc-1];
		
		return clusterings[k-1]; 
	}
	/**
	 * @param emResultList
	 */
	public void setClusterings(List list) {
		if (list==null) return;
		nmqc   = list.size();
		clusterings = new EMResult[nmqc];
		list.toArray(clusterings);
	}
	/**
	 * @param emResults
	 */
	public void setClusterings(EMResult[] results) {
		if (results==null) return;
		nmqc   = results.length;
		clusterings = results;
	}
	/**
	 * Set the indices of the Qtl which are implied in
	 * this result.
	 * 
	 * @param qtlIdx the indices of the Qtl.
	 */
	public void setQtlIdx(int[] qtlIdx) {
		this.qtlIdx = qtlIdx;
	}
	/**
	 * Set the observed positions of the Qtl.
	 * @param xx
	 */
	public void setX(Double[] xx) {
		if (xx == null) return;
		
		for(int i=0;i<nqtl;i++)
			x[i] = xx[i].doubleValue(); 
	}
	/**
	 * Set the standard deviations of the Qtl.
	 * @param xx
	 */
	public void setSD(Double[] sd) {
		if (sd == null) return;
		
		for(int i=0;i<nqtl;i++)
			this.sd[i] = sd[i].doubleValue(); 
	}
	/**
	 * 
	 * @param tree
	 */
	public void setTree(Tree tree) {
		this.tree = tree;
	}
	/**
	 * This methods returns the best models for each criterion
	 * defined to select the optimal number of clusters per
	 * trait.
	 * 
	 * @return
	 */
	public EMResult getBestClustering(int criterion) {
		return EMCriteria.getBestResult(clusterings, criterion);
	}
	/**
	 * Returns the name of the trait cluster for this result.
	 * @return the name of the trait cluster.
	 */
	public String getTraitName() {
		return trait;
	}
	/**
	 * @return
	 */
	public int[] getKValues() {
		if (clusterings==null) return null;
		/*
		* Count the number of models for which
		* the loglikelihood increases.
		*/
		int nk=0;
		for(int i=0;i<clusterings.length;i++)
			if (clusterings[i]!=null) nk++;
		int[] K = new int[nk];
		for(int j=0,i=0;i<clusterings.length;i++)
			if (clusterings[i]!=null)
				K[j++]=clusterings[i].k;
		return K;
	}
	/**
	 * @return
	 */
	public Iterator getCriteria() {
		return EMCriteria.getCriteria();
	}
	/**
	 * @param criterion
	 * @return
	 */
	public double[] getCritValues(String criterion, int[] K) {
		if (clusterings==null) return null;
		int cidx = EMCriteria.getCriterionIdx(criterion);
		double[] values = new double[K.length];
		for(int i=0;i<clusterings.length;i++) {
			if (clusterings[i]!=null) {
				for(int j=0;j<K.length;j++) {
					if (K[j]==clusterings[i].k) {
						values[j]=clusterings[i].getCriterion(cidx);
						break;
					}
				}	
			}
		}
		return values;
	}
	/**
	 * @param criterion
	 * @return
	 */
	public double[] getCritDelta(String criterion, int[] K) {
		if (clusterings==null) return null;
		int cidx = EMCriteria.getCriterionIdx(criterion);
		double[] delta = new double[K.length];
		double best    = getBestClustering(cidx).getCriterion(cidx);
		for(int i=0;i<clusterings.length;i++) {
			if (clusterings[i]!=null) {
				for(int j=0;j<K.length;j++) {
					if (K[j]==clusterings[i].k) {
						delta[j]=clusterings[i].getCriterion(cidx)-best;
					}
				}
			}
		}
		return delta;
	}
	/**
	 * 
	 * @param criterion
	 * @return
	 */
	public double[] getCritWeights(String criterion, int[] K) {
		if (clusterings==null) return null;
		int cidx = EMCriteria.getCriterionIdx(criterion);
		double[] weights = new double[K.length];
		double   best    = getBestClustering(cidx).getCriterion(cidx);
		double   sum     = 0.0;
		for(int i=0;i<clusterings.length;i++) {
			if (clusterings[i]!=null) {
				for(int j=0;j<K.length;j++) {
					if (K[j]==clusterings[i].k) {
						weights[j]=Math.exp(-0.5*(clusterings[i].getCriterion(cidx)-best));
						sum+=weights[j];
					}
				}
			}
		}
		for(int i=0;i<K.length;i++)weights[i] /= sum;
		return weights;
	}
	/**
	 * @param best
	 * @return
	 */
	public double[] getUSD(EMResult best, int kmin, int kmax, String criterion) {
		int[] K      = getKValues();
		int[] nsd    = new int[best.k];
		double[] w   = getCritWeights(criterion, K);
		double[] usd = new double[best.k];
		double[] xsd = new double[best.n];
		double[] wxp = new double[best.n];
		// First compute the model averaging predicted values.
		for(int j=0;j<K.length;j++) {
			if (K[j] >= kmin && K[j] <= kmax) {
				EMResult result = this.getClustering(K[j]);
				double[] xp     = result.getXBestPred();
				for(int i=0;i<result.n;i++)
					wxp[i] += w[j]*xp[i];
			}
		}
		// Then compute the model averaging standard deviation of the predicted values
		for(int j=0;j<K.length;j++) {
			if (K[j] >= kmin && K[j] <= kmax) {
				EMResult result = this.getClustering(K[j]);
				for(int i=0;i<result.n;i++) {
					int    xk = result.getXPredIdx(i);
					double sd = (wxp[i]-result.getMu(xk));
					sd *= sd;
					sd += result.ocov[xk][xk];
					xsd[i] += w[j]*Math.sqrt(sd);
				}
			}
		}
		for(int i=0;i<best.n;i++) {
			int xk   = best.getXPredIdx(i);
			nsd[xk]++;
			usd[xk] += xsd[i]; 
		}
		for(int i=0;i<best.k;i++) usd[i] /= (double)nsd[i];
		
		return usd;
	}
	/**
	 * @return
	 */
	public Tree getQtlTree() {
		return tree;
	}
	/**
	 * @return
	 */
	public Qtl[] getMetaQtl(int k) {
		EMResult result = getClustering(k);
		if (result == null) return null;
		Qtl[] qtls      = QtlAdapter.adapt(result);
		if (qtls != null) {
			// Set the trait
			for(int i=0;i<qtls.length;i++) {
				qtls[i].name       = new String("MQTL_"+(i+1)+"_"+trait);
				qtls[i].trait 	   = qtls[i].new QTLTrait();
				qtls[i].trait.name = trait;
			}
		}
		return qtls;
	}
}
