package dr.evomodel.substmodel;

import dr.inference.model.Variable;
import dr.inference.model.Parameter;
import dr.inference.model.Bounds;

/**
 * @author Chieh-Hsi Wu
 *
 * Model averaging of nucleotide substitution model
 *
 */
public class NtdBMA extends AbstractNucleotideModel{

    private Variable<Double> logKappa = null;
    private Variable<Double> logTN = null;
    private Variable<Double> logAC = null;
    private Variable<Double> logAT = null;
    private Variable<Double> logGC = null;
    private Variable<Double> logGT = null;
    private Variable<Integer> modelChoose = null;

    public NtdBMA(
            Variable<Double> logKappa,
            Variable<Double> logTN,
            Variable<Double> logAC,
            Variable<Double> logAT,
            Variable<Double> logGC,
            Variable<Double> logGT,
            Variable<Integer> modelChoose,
            FrequencyModel freqModel){
        super("NucleotideBMA", freqModel);


        addVariable(logKappa);
        logKappa.addBounds(new Parameter.DefaultBounds(Double.POSITIVE_INFINITY, Double.NEGATIVE_INFINITY, 1));
        this.logKappa = logKappa;

        addVariable(logTN);
        logTN.addBounds(new Parameter.DefaultBounds(Double.POSITIVE_INFINITY, Double.NEGATIVE_INFINITY, 1));
        this.logTN = logTN;

        addVariable(logAC);
        logAC.addBounds(new Parameter.DefaultBounds(Double.POSITIVE_INFINITY, Double.NEGATIVE_INFINITY, 1));
        this.logAC = logAC;

        addVariable(logAT);
        logAT.addBounds(new Parameter.DefaultBounds(Double.POSITIVE_INFINITY, Double.NEGATIVE_INFINITY, 1));
        this.logAT = logAT;

        addVariable(logGC);
        logGC.addBounds(new Parameter.DefaultBounds(Double.POSITIVE_INFINITY, Double.NEGATIVE_INFINITY, 1));
        this.logGC = logGC;

        addVariable(logGT);
        logGT.addBounds(new Parameter.DefaultBounds(Double.POSITIVE_INFINITY, Double.NEGATIVE_INFINITY, 1));
        this.logGT = logGT;

        addVariable(modelChoose);
        modelChoose.addBounds(new Bounds.Int(modelChoose, 0, 1));

        
        this.modelChoose = modelChoose;
    }

    protected void setupRelativeRates() {

        //rate AG value
        relativeRates[1] = Math.exp(logKappa.getValue(0)+modelChoose.getValue(0)*logTN.getValue(0));
        //rate CT value
        relativeRates[4] = Math.exp(logKappa.getValue(0)-modelChoose.getValue(0)*logTN.getValue(0));
        //rate AC value
        relativeRates[0] = Math.exp(modelChoose.getValue(1)*logAC.getValue(0));
        //rate AT value
        relativeRates[2] = Math.exp(modelChoose.getValue(1)*logAT.getValue(0));
        //rate GC value
        relativeRates[3] = Math.exp(modelChoose.getValue(1)*logGC.getValue(0));
        //rate GT value
        relativeRates[5] = Math.exp(modelChoose.getValue(1)*logGT.getValue(0));

    }

}

@@ -2919,11 +2919,1180 @@
    %7D%0D%0A%0D%0A
+    protected void handleVariableChangedEvent(Variable variable, int index, Parameter.ChangeType type) %7B%0D%0A        boolean changeRateMatrix = false;%0D%0A        %0D%0A        if(variable == modelChoose)%7B%0D%0A            //changing the substitution model%0D%0A            changeRateMatrix = true;%0D%0A        %7Delse if(variable == logKappa)%7B%0D%0A            //changing kappa value (which is in all models)%0D%0A            changeRateMatrix = true;%0D%0A        %7Delse if(variable == logTN && modelChoose.getValue(0) == 1)%7B%0D%0A            //changing the tn value when the current model is TN.%0D%0A            changeRateMatrix = true;%0D%0A        %7Delse if(variable == logAC %7C%7C%0D%0A                variable == logAT %7C%7C%0D%0A                variable == logGC %7C%7C%0D%0A                variable == logGT &&%0D%0A                modelChoose.getValue(1) == 1 && modelChoose.getValue(0) == 1)%7B%0D%0A            //Changing any one of A%3C-%3EC, A%3C-%3ET, G%3C-%3EC, G%3C-%3ET rates while the current model is GTR.%0D%0A            changeRateMatrix = true;%0D%0A        %7D%0D%0A        // relativeRates changed%0D%0A        if(changeRateMatrix)%7B%0D%0A            super.handleVariableChangedEvent(variable, index, type);%0D%0A        %7D%0D%0A        changeRateMatrix = false;%0D%0A    %7D%0D%0A%0D%0A
 %7D%0D%0A

Handle the variable change more efficiently.