C++ Genome::getBestRecombinantTriplet方法代码示例

void FullRun::perform() {
    Run::perform(); // process common data

    verifyData();
    setup();
    loadPTable(pTableFile);


    Interface::instance() << Interface::SEPARATOR << endl;
    Interface::instance().showLog(true);

    progress();

    Interface::instance() << Interface::SEPARATOR << endl;
    Interface::instance().showLog(true);

    if (!cloAllBpCalculated && !cloNoBpCalculated && !cloNo3sRecFile) {
        /* We calculate breakpoints for the best triplets */
        for (unsigned long i = 0; i < childDataset->getSize(); i++) {
            Genome* child = childDataset->getGenome(i);
            Triplet* bestRecTriplet = child->getBestRecombinantTriplet();
            if (bestRecTriplet != NULL) {
                recordRecombinantTriplet(bestRecTriplet);
            }
        }
    }


    if (isRandomMode && randomLoopNum > 1) {
        Interface::instance() << "Recombination sensitivity: "
                << recombinationSensitivity << "/" << randomLoopNum << "\n";
        Interface::instance().showOutput(true);
        
    } else {
        displayResult();
        Interface::instance() << Interface::SEPARATOR << endl;
        Interface::instance().showLog(true);
        savePValHistogram();
    }


    /* Close all files */
    if (fileSkippedTriplets)
        fileSkippedTriplets->close();
    if (fileRecombinants)
        fileRecombinants->close();
    if (fileLongRecombinants != NULL)
        fileLongRecombinants->close();
}

//.........这里部分代码省略.........
        Genome* mum = NULL;
        Triplet* triplet = NULL;

        for (unsigned long childIndex = 0; !isStoped && childIndex < childDataset->getSize(); childIndex++) {
            child = childDataset->getGenome(childIndex);
            if (child->isActive()) {
                activeChildCounter++;
            } else {
                continue;
            }

            activeDadCounter = 0;
            for (unsigned long dadIndex = 0; !isStoped && dadIndex < parentDataset->getSize(); dadIndex++) {
                dad = parentDataset->getGenome(dadIndex);
                if (dad->isActive()) {
                    activeDadCounter++;
                }
                if (!dad->isActive() || dad == child) {
                    continue;
                }

                /* Show progress counter */
                double currentLoop = static_cast<double> (activeChildCounter - 1)
                        * static_cast<double> (activeParentNum)
                        + static_cast<double> (activeDadCounter - 1);
                showProgress(currentLoop, false);

                for (unsigned long mumIndex = 0; !isStoped && mumIndex < parentDataset->getSize(); mumIndex++) {
                    mum = parentDataset->getGenome(mumIndex);
                    if (!mum->isActive() || mum == child || mum == dad) continue;

                    triplet = new Triplet(dad, mum, child);

                    /* Let's see if we can calculate the exact P-value */
                    bool hasPValue = triplet->calculatePVal(cloUseSiegmundApprox);
                    if (!hasPValue) {
                        if (fileSkippedTriplets) {
                            /* No need to open file, writeLine() will do that automatically */
                            fileSkippedTriplets->writeLine(triplet->toString());
                        }
                        delete triplet;
                        continue;
                    }
                    numSkipped -= 1.0;

                    if (triplet->isPValApproximated()) {
                        numApproximated += 1.0;
                    } else {
                        numComputedExactly += 1.0;
                    }

                    double pValue = triplet->getPVal();
                    addPValIntoHistogram(pValue);

                    if (pValue < minPVal)
                        minPVal = pValue;

                    /* Let's see if the P-value is significant */
                    if (stats::correction::dunnSidak(pValue, getNumTripletsForCorrection()) < cloRejectThreshold) {
                        numRecombinantTriplets += 1.0;
                        if (child->getRecombinantType() == Genome::Na_REC) {
                            child->setRecombinantType(Genome::SHORT_REC);
                        }

                        /* Go into this loop either if we're doing ALL breakpoint or
                         * NO breakpoint if we are doing no breakpoint, then the 
                         * breakpoint calculations is skipped. */
                        if (cloAllBpCalculated || cloNoBpCalculated) {
                            recordRecombinantTriplet(triplet);
                        }

                        /* Now check if this is the *best* recombinant (meaning lowest p-value)
                         * and if it is, record it */
                        if (child->getBestRecombinantTriplet() == NULL
                                || pValue < child->getBestRecombinantTriplet()->getPVal()) {
                            child->setBestRecombinantTriplet(triplet);
                        } else {
                            delete triplet;
                        }

                        if (cloStopAtFirstRecombinant) {
                            /* Finish progressing early */
                            isStoped = true;
                        }

                    } else {
                        delete triplet;
                    }
                }
            }
        }

        /* Progressing finished */
        showProgress(0.0, true);

        if (numRecombinantTriplets > 0) {
            recombinationSensitivity++;
        }
    }
}