package normal2;
import java.io.File;
import java.io.IOException;
import java.io.PrintWriter;
import java.util.Arrays;
import org.apache.commons.math.random.RandomDataImpl;
import org.apache.commons.math.random.RandomGenerator;
import org.apache.commons.math.random.Well44497b;
import common.FragmentPartitionPrior;
import common.FragmentPenalty;
import common.PoissonTuningDistribution;
import common.PowerFragmentPenalty;
public class RunMCMC {
public static void main(String[] args) throws IOException {
int n = 20;
String dataFile = "/home/adu/Dropbox/research/vardeman/mc-on-partition/data/jdata1/data-std.bin";
String iapFile = "/home/adu/Dropbox/research/vardeman/mc-on-partition/data/jdata1/apt.bin";
double[] data = dclong.io.BinaryReader.readDouble(dataFile, 0, 0, n);
int[] initArrayPartition = dclong.io.BinaryReader.readInt(iapFile,0,0,n);
// String outputFile = "/home/adu/Dropbox/research/vardeman/mc-on-partition/data/jdata2/r-0.5-6.0-1.txt";
// if(new File(outputFile).exists()){
// System.err.println("The specified output file already exists. Please use a new one.");
// System.exit(1);
// }
double power = 0.5;
double powerCoefficient = 6;
double poisson = 0.2;
double alpha = 2;
double beta = 0.2;
int warmupStep = 10000;
int numberOfDraws = 100000;
String output = "PowerFragmentPenalty Settings:\npower="+power+", coefficient="+powerCoefficient+"\n";
output += "\nVariance Prior Settings:\nalpha="+alpha+", beta="+beta+"\n";
output += "\nTuning Settings:\nG: uniform, F: 1+Poisson("+poisson+")\n";
output += "\nRunning Settings:\n";
output += "Initial partition in array representation:\n"+Arrays.toString(initArrayPartition)+"\n";
output += "Initial partition in list representation:\n"+dclong.util.Arrays.splitArray(initArrayPartition)+"\n";
output += "burn-in="+warmupStep+", draws="+numberOfDraws+"\n";
//define partition prior
// FragmentPenalty fragmentPenalty = new PowerFragmentPenalty(0.5,3);//seems good
FragmentPenalty fragmentPenalty = new PowerFragmentPenalty(power,powerCoefficient);
//define partition prior
FragmentPartitionPrior fragPrior = new FragmentPartitionPrior(fragmentPenalty);
//define tuning distribution
RandomGenerator rg = new Well44497b();
RandomDataImpl rng = new RandomDataImpl(rg);
common.TuningDistribution tuningDistribution = new PoissonTuningDistribution(rng,poisson);
//initial jumping weights
double[] initJumpWeight = new double[n];
java.util.Arrays.fill(initJumpWeight, 10);
//------------------------------------------------------------------------------------
MCMC mcmc = new MCMC(rng, data, initArrayPartition,
initJumpWeight,fragPrior,tuningDistribution,alpha,beta);
// System.out.println("1 class complexity is: "+mcmc.classComplexity());
dclong.util.Timer timer = new dclong.util.Timer();
timer.start();
mcmc.warmup(warmupStep,false);
timer.stop();
output += "\nTime used for warmup is: " + timer.seconds()+" seconds.\n";
output += "Warmup jumping ratio: " + mcmc.getJumpRatio() + "\n";
// timer.start();
// mcmc.warmup(10, false);
// timer.stop();
// timer.printSeconds("warmup without auto adjusting");
// System.out.println("Jump ratio: "+mcmc.getJumpRatio());
// System.out.println(dclong.util.Arrays.sum(mcmc.getJumpWeights())-100*20);
// dclong.io.Print.print(dclong.util.Arrays.cbind(dclong.util.Arrays.sequence(0, 19, 1),mcmc.getJumpWeights()));
// dclong.io.Print.print(mcmc.getProbability());
//
timer.start();
mcmc.run(numberOfDraws);
timer.stop();
output += "\nTime used for simulation is: " + timer.seconds() + " seconds.\n";
output += "Simulation jumping ratio: " + mcmc.getJumpRatio() + "\n";
//write partitions into disk
// mcmc.writeFP("/home/adu/Dropbox/research/vardeman/mc-on-partition/simulation/draws-ap1-f13-s0.1.bin");
output += "Unshrinked number of partitions: " + mcmc.getFps().size()+"\n";
// timer.start();
// common.FreqPartition center = common.FreqPartition.center(mcmc.getFps(), true);
// timer.end();
// output += "Time used for finding the center is: " + timer.seconds() + " seconds.\n";
// output += "Shrinked number of partitions: " + mcmc.getFps().size()+"\n";
// output += "\nThe center in array representation is:\n"+center+"\n";
// output += "The center in list representation is:\n"+dclong.util.Arrays.splitArray(center.getPartition())+"\n";
//
timer.start();
common.FreqPartition.reduce(mcmc.getFps());
timer.end();
output += "Time used for reduction is: " + timer.seconds() + " seconds.\n";
output += "Shrinked number of partitions: " + mcmc.getFps().size() + "\n";
common.FreqPartition mode = common.FreqPartition.mode(mcmc.getFps());
output += "\nThe mode in array representation is:\n"+mode+"\n";
output += "The mode in list representation is:\n"+dclong.util.Arrays.splitArray(mode.getPartition())+"\n";
System.out.println(output);
// PrintWriter pw = new PrintWriter(outputFile);
// pw.write(output);
// pw.close();
}
}