Java Randomizer.randomChoicePDF方法代码示例

import beast.util.Randomizer; //导入方法依赖的package包/类
/**
 * perform the actual sequence generation
 *
 * @return alignment containing randomly generated sequences for the nodes in the
 *         leaves of the tree
 */
public void simulate() {
    Node root = m_tree.getRoot();


    double[] categoryProbs = m_siteModel.getCategoryProportions(root);
    int[] category = new int[m_sequenceLength];
    for (int i = 0; i < m_sequenceLength; i++) {
        category[i] = Randomizer.randomChoicePDF(categoryProbs);
    }

    double[] frequencies = m_siteModel.getSubstitutionModel().getFrequencies();
    int[] seq = new int[m_sequenceLength];
    for (int i = 0; i < m_sequenceLength; i++) {
        seq[i] = Randomizer.randomChoicePDF(frequencies);
    }


    //alignment.setDataType(m_siteModel.getFrequencyModel().getDataType());

    traverse(root, seq, category);

}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
/**
 * recursively walk through the tree top down, and add sequence to alignment whenever
 * a leave node is reached.
 *
 * @param node           reference to the current node, for which we visit all children
 * @param parentSequence randomly generated sequence of the parent node
 * @param category       array of categories for each of the sites
 * @param alignment
 */
void traverse(Node node, int[] parentSequence, int[] category) {
    for (int childIndex = 0; childIndex < 2; childIndex++) {
        Node child = (childIndex == 0 ? node.getLeft() : node.getRight());
        for (int i = 0; i < m_categoryCount; i++) {
            getTransitionProbabilities(m_tree, child, i, m_probabilities[i]);
        }

        int[] seq = new int[m_sequenceLength];
        double[] cProb = new double[m_stateCount];
        for (int i = 0; i < m_sequenceLength; i++) {
            System.arraycopy(m_probabilities[category[i]], parentSequence[i] * m_stateCount, cProb, 0, m_stateCount);
            seq[i] = Randomizer.randomChoicePDF(cProb);
        }

        if (child.isLeaf()) {
            sequenceInput.setValue(intArray2Sequence(seq, child), this);
        } else {
            traverse(child, seq, category);
        }
    }
}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
/**
 * recursively walk through the tree top down, and add sequence to alignment whenever
 * a leave node is reached.
 *
 * @param node           reference to the current node, for which we visit all children
 * @param parentSequence randomly generated sequence of the parent node
 * @param category       array of categories for each of the sites
 * @param alignment
 * @
 */
void traverse(Node node, int[] parentSequence, int[] category, Alignment alignment)  {
    for (int childIndex = 0; childIndex < 2; childIndex++) {
        Node child = (childIndex == 0 ? node.getLeft() : node.getRight());
        for (int i = 0; i < m_categoryCount; i++) {
            getTransitionProbabilities(m_tree, child, i, m_probabilities[i]);
        }

        int[] seq = new int[m_sequenceLength];
        double[] cProb = new double[m_stateCount];
        for (int i = 0; i < m_sequenceLength; i++) {
            System.arraycopy(m_probabilities[category[i]], parentSequence[i] * m_stateCount, cProb, 0, m_stateCount);
            seq[i] = Randomizer.randomChoicePDF(cProb);
        }

        if (child.isLeaf()) {
            alignment.sequenceInput.setValue(intArray2Sequence(seq, child), alignment);
        } else {
            traverse(child, seq, category, alignment);
        }
    }
}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
protected int getDeathNode(String[] aliveNodes, int[] traits, int numberOfLangs) {
	double totalCognates = 0.0;
	for (int n : traits) {
		totalCognates += n;
	}
	double[] probs = new double[numberOfLangs];
	for (int i = 0; i < numberOfLangs; i++) {
		probs[i] = traits[i] / totalCognates;
	}
	return Randomizer.randomChoicePDF(probs);
}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
/**
 * perform the actual sequence generation
 *
 * @return alignment containing randomly generated sequences for the nodes in the
 *         leaves of the tree
 * @
 */
public Alignment simulate()  {
    Node root = m_tree.getRoot();


    double[] categoryProbs = m_siteModel.getCategoryProportions(root);
    int[] category = new int[m_sequenceLength];
    for (int i = 0; i < m_sequenceLength; i++) {
        category[i] = Randomizer.randomChoicePDF(categoryProbs);
    }

    double[] frequencies = m_siteModel.getSubstitutionModel().getFrequencies();
    int[] seq = new int[m_sequenceLength];
    for (int i = 0; i < m_sequenceLength; i++) {
        seq[i] = Randomizer.randomChoicePDF(frequencies);
    }


    Alignment alignment = new Alignment();
    alignment.userDataTypeInput.setValue(m_data.get().getDataType(), alignment);
    alignment.setID("SequenceSimulator");

    traverse(root, seq, category, alignment);


    return alignment;
}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
/**
 * Traverse a tree, simulating a sequence alignment down it.
 *
 * @param node Node of the tree
 * @param parentSequence Sequence at the parent node in the tree
 * @param categories Mapping from sites to categories
 * @param transitionProbs
 * @param regionAlignment
 */
private void traverse(Node node,
        int[] parentSequence,
        int[] categories, double[][] transitionProbs,
        int[][] regionAlignment) {

    for (Node child : node.getChildren()) {

        // Calculate transition probabilities
        for (int i=0; i<siteModel.getCategoryCount(); i++) {
            siteModel.getSubstitutionModel().getTransitionProbabilities(
                    child, node.getHeight(), child.getHeight(),
                    siteModel.getRateForCategory(i, child),
                    transitionProbs[i]);
        }

        // Draw characters on child sequence
        int[] childSequence = new int[parentSequence.length];
        int nStates = dataType.getStateCount();
        double[] charProb = new double[nStates];
        for (int i=0; i<childSequence.length; i++) {
            int category = categories[i];
            System.arraycopy(transitionProbs[category],
                    parentSequence[i]*nStates, charProb, 0, nStates);
            childSequence[i] = Randomizer.randomChoicePDF(charProb);
        }

        if (child.isLeaf()) {
            System.arraycopy(childSequence, 0,
                    regionAlignment[child.getNr()], 0, childSequence.length);
        } else {
            traverse(child, childSequence,
                    categories, transitionProbs,
                    regionAlignment);
        }
    }
}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
/**
 * Traverse a marginal tree simulating a region of the sequence alignment
 * down it.
 * 
 * @param node Node of the marginal tree
 * @param parentSequence Sequence at the parent node in the marginal tree
 * @param categories Mapping from sites to categories
 * @param transitionProbs
 * @param regionAlignment 
 */
private void traverse(Node node,
        int[] parentSequence,
        int[] categories, double[][] transitionProbs,
        int[][] regionAlignment) {
    
    for (Node child : node.getChildren()) {

        // Calculate transition probabilities
        for (int i=0; i<siteModel.getCategoryCount(); i++) {
            siteModel.getSubstitutionModel().getTransitionProbabilities(
                    child, node.getHeight(), child.getHeight(),
                    siteModel.getRateForCategory(i, child),
                    transitionProbs[i]);
        }
        
        // Draw characters on child sequence
        int[] childSequence = new int[parentSequence.length];
        int nStates = dataType.getStateCount();
        double[] charProb = new double[nStates];
        for (int i=0; i<childSequence.length; i++) {
            int category = categories[i];
            System.arraycopy(transitionProbs[category],
                    parentSequence[i]*nStates, charProb, 0, nStates);
            childSequence[i] = Randomizer.randomChoicePDF(charProb);
        }
        
        if (child.isLeaf()) {
            System.arraycopy(childSequence, 0,
                    regionAlignment[child.getNr()], 0, childSequence.length);
        } else {
            traverse(child, childSequence,
                    categories, transitionProbs,
                    regionAlignment);
        }
    }
}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
/**
 * Select affected region when region endpoint restriction is in place.
 *
 * @param conv Conversion object where these sites are stored.
 * @return log probability density of chosen attachment.
 */
protected double drawAffectedRegionRestricted(Conversion conv) {
    int locusIdx = Randomizer.randomChoicePDF(relativeLocusSizes);
    Locus locus = acg.getLoci().get(locusIdx);

    conv.setLocus(locus);
    conv.setStartSite(0);
    conv.setEndSite(locus.getSiteCount()-1);

    return Math.log(relativeLocusSizes[locusIdx]);
}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
public Tree mutateOverTreeBorrowing(Tree base) throws Exception {
	Double[] events = getEvents(base);
	setSubTreeLanguages(base.getRoot(), (Sequence) base.getRoot().getMetaData("lang"));
	// Get root node.
	ArrayList<Node> aliveNodes = new ArrayList<Node>();
	String[] stringAliveNodes = {};
	int[] traits = {};
	int numberOfLangs = 0;
	Double totalRate = null;
	int idx, ind;
	double[] probs = new double[3];
	for (int i = 0; i < events.length - 1; i++) {
		if (events[i] == 0.0) {
			break;
		}
		aliveNodes = getAliveNodes(base, events[i+1]);
		stringAliveNodes = getSequences(aliveNodes);
		numberOfLangs = stringAliveNodes.length;
		traits = getBirths(stringAliveNodes, numberOfLangs);
		totalRate = totalRate(stringAliveNodes, traits, numberOfLangs);
		Double t = events[i] - Randomizer.nextExponential(totalRate);
		//System.out.println();
		//System.out.println("On branch event: " + (i+1)+ " out of " + (events.length/2) + ". Next event at " + events[i+1]);
		while (t > events[i+1]) {
			//System.out.print("\r"+t);
			probs = BorrowingProbs(stringAliveNodes, totalRate,traits, numberOfLangs);
			Integer choice = Randomizer.randomChoicePDF(probs);
			// Birth.
			if (choice == 0) {
				idx = Randomizer.nextInt(stringAliveNodes.length);
				stringAliveNodes[idx] = stringAliveNodes[idx] + "1";
				stringAliveNodes = addEmptyTrait(stringAliveNodes, idx);
				traits[idx]++;
				// Death.
			} else if (choice == 1) {
				idx = getDeathNode(stringAliveNodes, traits, numberOfLangs);
				// Find random alive trait, and kill it.
				ind = getRandomBirthIndex(stringAliveNodes[idx], traits[idx]);
				if (noEmptyTraitCheck(stringAliveNodes[idx]) && ind > -1) {
					stringAliveNodes[idx] = replaceCharAt(stringAliveNodes[idx], ind, Integer.toString(0));
					traits[idx]--;
				}
				// Borrowing.
			} else if (choice == 2) {
				if (aliveNodes.size() > 1) {
					// Pick two distinct languages at random.
					int[] bN = getBorrowingNodes(stringAliveNodes, traits, numberOfLangs);

					if (localDist(aliveNodes.get(bN[0]), aliveNodes.get(bN[1]))) {
						// Randomly iterate through language and find a 1.
						ind = getRandomBirthIndex(stringAliveNodes[bN[0]], traits[bN[0]]);
						// If recieving language is going 0 -> 1.
						if (ind > 1 && stringAliveNodes[bN[1]].charAt(ind) == '0') {
							traits[bN[1]]++;
						}
						// Give the 1 to the receiving language.
						stringAliveNodes[bN[1]] = replaceCharAt(stringAliveNodes[bN[1]], ind, Integer.toString(1));
					}
				}
			}
			totalRate = totalRate(stringAliveNodes, traits, numberOfLangs);
			t -= Randomizer.nextExponential(totalRate);
		}
		setLangs(aliveNodes, stringAliveNodes);
	}
	return base;
}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
public Tree mutateOverTreeBorrowing(Tree base) throws Exception {
	Double[] events = getEvents(base);
	setSubTreeLanguages(base.getRoot(), (Sequence) base.getRoot().getMetaData("lang"));
	// Get root node.
	ArrayList<Node> aliveNodes = new ArrayList<Node>();
	String[] stringAliveNodes = {};
	int[] traits = {};
	int numberOfLangs = 0;
	// Get first event.
	Double totalRate = null;
	//Double t = treeHeight - Randomizer.nextExponential(totalRate);
	// Variable declarations.
	int idx;
	double[] probs;
	for (int i = 0; i < events.length - 1; i++) {
		if (events[i] == 0.0) {
			break;
		}
		aliveNodes = getAliveNodes(base, events[i+1]);
		stringAliveNodes = getSequences(aliveNodes);
		numberOfLangs = stringAliveNodes.length;
		traits = getBirths(stringAliveNodes, numberOfLangs);
		totalRate = totalRate(stringAliveNodes, traits, numberOfLangs);
		Double t = events[i] - Randomizer.nextExponential(totalRate);
		//System.out.println();
		//System.out.println("On branch event: " + (i+1)+ " out of " + (events.length/2) + ". Next event at " + events[i+1]);
		while (t > events[i+1]) {
			//System.out.print("\r"+t);
			// Return array of event probabilities and pick one.
			probs = BorrowingProbs(stringAliveNodes, totalRate, traits, numberOfLangs);
			Integer choice = Randomizer.randomChoicePDF(probs);
			// Mutate.
			if (choice == 0) {
				// Pick a random node at time t.
				idx = Randomizer.nextInt(stringAliveNodes.length);
				// Pick a random position in language.
				int pos = Randomizer.nextInt(stringAliveNodes[idx].length());
				int currentTrait = Character.getNumericValue(stringAliveNodes[idx].charAt(pos));
				// If death and noEmptyTraitCheck fails.
				if (currentTrait == 1 && (!noEmptyTraitCheck(stringAliveNodes[idx]))) {
					stringAliveNodes[idx] = stringAliveNodes[idx];
				} else {
					stringAliveNodes[idx] = replaceCharAt(stringAliveNodes[idx], pos, Integer.toString((1 - currentTrait)));
				} 
				// Change trait structure
				if (stringAliveNodes[idx].charAt(pos) == '1') {
					traits[idx]++;
				} else {
					traits[idx]--;
				}

				// Borrow.
			} else if (choice == 1) {
				if (aliveNodes.size() > 1) {
					// Pick two distinct languages at random.
					int[] bN = getBorrowingNodes(stringAliveNodes, traits, numberOfLangs);

					if (localDist(aliveNodes.get(bN[0]), aliveNodes.get(bN[1]))) {
						// Randomly iterate through language and find a 1.
						int ind = getRandomBirthIndex(stringAliveNodes[bN[0]], traits[bN[0]]);
						// If recieving language is going 0 -> 1.
						if (ind > -1 && stringAliveNodes[bN[1]].charAt(ind) == '0') {
							traits[bN[1]]++;
						}
						// Give the 1 to the receiving language.
						stringAliveNodes[bN[1]] = replaceCharAt(stringAliveNodes[bN[1]], ind, Integer.toString(1));
					}
				}
			}
			totalRate = totalRate(stringAliveNodes, traits, numberOfLangs);
			t -= Randomizer.nextExponential(totalRate);
		}
		setLangs(aliveNodes, stringAliveNodes);
	}
	return base;
}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
public double proposal(){
    //Get the pointer and the list of unique values
    DPPointer pointers = pointersInput.get(this);
    ParameterList paramList = xListInput.get(this);

    //Randomly pick an index to update, gets it's current value and its position in the parameter list
    int dimPointer = pointers.getDimension();
    int index = Randomizer.nextInt(dimPointer);
    RealParameter currVal = paramList.getParameter(pointers.indexInList(index,paramList));
    int listIndex = paramList.indexOf(currVal);


    //Distribution lik = likelihoodInput.get();

    //Get the dimension of the parameter
    int dimValue = paramList.getParameterDimension();

    //Count the number of items in each cluster but excluding the one about to be updated
    int[] clusterCounts = dpVal.getClusterCounts();
    clusterCounts[listIndex] =  clusterCounts[listIndex]-1;

    try{

        //Generate a sample of proposals
        QuietRealParameter[] preliminaryProposals = sampleFromBaseDistribution(dimValue);
        int dimList = paramList.getDimension();
        int i;
        double concVal =dp.getConcParameter();
        
        double[] fullConditional = new double[dimList+sampleSize];

        for(i = 0; i < dimList; i++){
            //n_{-index,i}/(n - 1 + alpha)
            fullConditional[i] = clusterCounts[i]/(dimPointer - 1 + concVal);
            //System.err.println("fullConditional[i]: "+fullConditional[i]+", val: "+paramList.getParameter(i));

        }
        for(; i < fullConditional.length; i++){
            //alpha/m/(n - 1 + alpha)
            fullConditional[i] = concVal/sampleSize/(dimPointer - 1 + concVal);
            //System.err.println("fullConditional[i]: "+fullConditional[i]+", val: "+i);
        }

        int proposedIndex = Randomizer.randomChoicePDF(fullConditional);
        //System.err.println("proposedIndex: "+proposedIndex);
        if(proposedIndex < dimList){
            //take up an existing value
            pointers.point(index, paramList.getParameter(proposedIndex));
        }else{
            //take up a new value
            pointers.point(index, preliminaryProposals[proposedIndex-dimList]);
            paramList.addParameter(preliminaryProposals[proposedIndex-dimList]);
        }

        /*for(i = 0; i < clusterCounts.length;i++){
            System.err.println(paramList.getParameter(i).getValue()+", cluster counts: "+clusterCounts[i]);
        }*/

        //If any cluster has no member then it is removed.
        for(i = 0; i < clusterCounts.length;i++){
            if(clusterCounts[i] ==0){
                paramList.removeParameter(i);
                break;
            }

        }

    }catch(Exception e){
            throw new RuntimeException(e);
        }
    return Double.POSITIVE_INFINITY;
}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
/**
 * Perform actual sequence simulation.
 */
private void simulate() {
    int nTaxa = tree.getLeafNodeCount();

    double[] categoryProbs = siteModel.getCategoryProportions(tree.getRoot());

    int nCategories = siteModel.getCategoryCount();
    int nStates = dataType.getStateCount();
    double[][] transitionProbs = new double[nCategories][nStates*nStates];

    int[][] alignment = new int[nTaxa][seqLength];

    int[] categories = new int[seqLength];
    for (int i=0; i<seqLength; i++)
        categories[i] = Randomizer.randomChoicePDF(categoryProbs);

    Node root = tree.getRoot();

    int[] parentSequence = new int[seqLength];
    double[] frequencies = siteModel.getSubstitutionModel().getFrequencies();
    for (int i=0; i<parentSequence.length; i++)
        parentSequence[i] = Randomizer.randomChoicePDF(frequencies);

    ancestralSeqStr = dataType.state2string(parentSequence);

    traverse(root, parentSequence,
            categories, transitionProbs,
            alignment);

    for (int leafIdx=0; leafIdx<nTaxa; leafIdx++) {
        String seqString = dataType.state2string(alignment[leafIdx]);

        String taxonName;
        if (tree.getNode(leafIdx).getID() != null)
            taxonName = tree.getNode(leafIdx).getID();
        else
            taxonName = "t" + leafIdx;

        sequenceInput.setValue(new Sequence(taxonName, seqString), this);
    }
}
 

import beast.util.Randomizer; //导入方法依赖的package包/类
/**
 * Perform actual sequence simulation.
 */
private void simulate(Locus locus) {
    Node cfRoot = acg.getRoot();
    int nTaxa = acg.getLeafNodeCount();
    
    double[] categoryProbs = siteModel.getCategoryProportions(cfRoot);

    int nCategories = siteModel.getCategoryCount();
    int nStates = dataType.getStateCount();
    double[][] transitionProbs = new double[nCategories][nStates*nStates];
    
    int[][] alignment = new int[nTaxa][locus.getSiteCount()];
    
    for (Region region : acg.getRegions(locus)) {
        int thisLength = region.getRegionLength();

        MarginalTree marginalTree = new MarginalTree(acg, region);
        
        int[] categories = new int[thisLength];
        for (int i=0; i<thisLength; i++)
            categories[i] = Randomizer.randomChoicePDF(categoryProbs);
        
        int[][] regionAlignment = new int[nTaxa][thisLength];
        
        Node thisRoot = marginalTree.getRoot();
        
        int[] parentSequence = new int[region.getRegionLength()];
        double[] frequencies = siteModel.getSubstitutionModel().getFrequencies();
        for (int i=0; i<parentSequence.length; i++)
            parentSequence[i] = Randomizer.randomChoicePDF(frequencies);

        traverse(thisRoot, parentSequence,
                categories, transitionProbs,
                regionAlignment);
        
        // DEBUG: Count differences
        int segsites = 0;
        for (int s=0; s<regionAlignment[0].length; s++) {
            int state = regionAlignment[0][s];
            for (int l=1; l<nTaxa; l++) {
                if (state != regionAlignment[l][s]) {
                    segsites += 1;
                    break;
                }
            }
        }
        System.out.println(segsites + " segregating sites in region " + region);
        
        copyToAlignment(alignment, regionAlignment, region);
    }
    
    for (int leafIdx=0; leafIdx<nTaxa; leafIdx++) {
        String sSeq = dataType.state2string(alignment[leafIdx]);
        String sTaxon = acg.getNode(leafIdx).getID();
        sequenceInput.setValue(new Sequence(sTaxon, sSeq), this);
    }
}