use of dr.stats.Variate in project beast-mcmc by beast-dev.
the class SkylineReconstructor method main.
public static void main(String[] argv) {
Variate x = null;
List<Variate> plots = new ArrayList<Variate>();
for (int i = 1; i <= 200; i++) {
String stem = "sim" + (i < 10 ? "00" : (i < 100 ? "0" : "")) + i;
try {
SkylineReconstructor skyline = new SkylineReconstructor(new File(stem + ".log"), new File(stem + ".trees"), 1000000, 200, 0.0, 150000, 0.0);
if (x == null) {
x = skyline.getXData();
}
plots.add(skyline.getYDataMean());
} catch (IOException e) {
e.printStackTrace();
} catch (ImportException e) {
e.printStackTrace();
} catch (TraceException e) {
e.printStackTrace();
}
if (i % 10 == 0) {
System.err.println("Read " + i);
}
}
for (int i = 0; i < x.getCount(); i++) {
System.out.print(x.get(i));
for (Variate y : plots) {
System.out.print("\t" + y.get(i));
}
System.out.println();
}
}
use of dr.stats.Variate in project beast-mcmc by beast-dev.
the class LineageCountThroughTime method getLTT.
public static Variate[] getLTT(String treeFile, double minTime, double maxTime, int binCount, // number of trees to skip
int skip) throws IOException, ImportException {
double delta = (maxTime - minTime) / (binCount - 1);
BufferedReader reader = new BufferedReader(new FileReader(treeFile));
String line = reader.readLine();
TreeImporter importer;
if (line.toUpperCase().startsWith("#NEXUS")) {
importer = new NexusImporter(reader);
} else {
importer = new NewickImporter(reader, false);
}
int state = 0;
while (importer.hasTree() && state < skip) {
importer.importNextTree();
state += 1;
}
// the age of the end of this group
List<double[]> branchingTimes = new ArrayList<double[]>();
state = 0;
while (importer.hasTree()) {
RootedTree tree = (RootedTree) importer.importNextTree();
double[] bt = new double[tree.getInternalNodes().size()];
int i = 0;
for (Node node : tree.getInternalNodes()) {
bt[i] = tree.getHeight(node);
i++;
}
Arrays.sort(bt);
branchingTimes.add(bt);
state += 1;
}
Variate[] bins = new Variate[binCount];
double height = 0;
double n = branchingTimes.get(0).length;
for (int k = 0; k < binCount; k++) {
bins[k] = new Variate.D();
if (height >= 0.0 && height <= maxTime) {
for (state = 0; state < branchingTimes.size(); state++) {
int index = 0;
while (index < branchingTimes.get(state).length && branchingTimes.get(state)[index] < height) {
index += 1;
}
double lineageCount = 1;
if (index < branchingTimes.get(state).length) {
lineageCount = n - index + 1;
}
bins[k].add(lineageCount);
}
}
height += delta;
}
Variate xData = new Variate.D();
Variate yDataMean = new Variate.D();
Variate yDataMedian = new Variate.D();
Variate yDataUpper = new Variate.D();
Variate yDataLower = new Variate.D();
double t = minTime;
for (Variate bin : bins) {
xData.add(t);
if (bin.getCount() > 0) {
yDataMean.add(bin.getMean());
yDataMedian.add(bin.getQuantile(0.5));
yDataLower.add(bin.getQuantile(0.025));
yDataUpper.add(bin.getQuantile(0.975));
} else {
yDataMean.add(Double.NaN);
yDataMedian.add(Double.NaN);
yDataLower.add(Double.NaN);
yDataUpper.add(Double.NaN);
}
t += delta;
}
return new Variate[] { xData, yDataMean };
}
use of dr.stats.Variate in project beast-mcmc by beast-dev.
the class KernelDensityEstimate method calculateDensity.
protected void calculateDensity(Variate data, int minimumBinCount) {
FrequencyDistribution frequency = calculateFrequencies(data, minimumBinCount);
xCoordinates = new Variate.D();
yCoordinates = new Variate.D();
double x = frequency.getLowerBound() - (frequency.getBinSize() / 2.0);
int extraEdgeCount = 0;
while (kde.pdf(x) > minDensity && x > lowerBoundary) {
x -= frequency.getBinSize();
extraEdgeCount += 1;
}
xCoordinates.add(x);
yCoordinates.add(0.0);
x += frequency.getBinSize();
int count = 0;
while (count < (frequency.getBinCount() + extraEdgeCount)) {
// ||
// (kde.pdf(x) > minDensity && x < upperBoundary)) {
xCoordinates.add(x);
yCoordinates.add(kde.pdf(x));
x += frequency.getBinSize();
count++;
}
// System.err.println("kde = " + kde.pdf(x));
while (kde.pdf(x) > minDensity) {
// System.err.println("add bit on end!!!");
xCoordinates.add(x);
yCoordinates.add(kde.pdf(x));
x += frequency.getBinSize();
}
xCoordinates.add(x);
yCoordinates.add(0.0);
//
//
// int extraBinsOnEdges = 5;
// double x = frequency.getLowerBound() - extraBinsOnEdges * frequency.getBinSize();
// for (int i = 0; i < frequency.getBinCount() + 2 * extraBinsOnEdges; i++) {
// double xMidPoint = x + (frequency.getBinSize() / 2.0);
// xData.add(xMidPoint);
// yData.add(kde.pdf(xMidPoint));
// x += frequency.getBinSize();
// }
}
Aggregations