Method calculateRatioUsingCG

src/GenomeCopyNumber.cpp:568–802

}

int GenomeCopyNumber::calculateRatioUsingCG (bool intercept, float minExpectedGC, float maxExpectedGC) { //returns 1 if the number of iteration is less than max; otherwise 0

     //try degree 3 and 4, SELECT THE ONE WITH LESS ITTERATIONS:
    int minDegreeToTest = 3;
    int maxDegreeToTest = 4;
    int selectedDegree=minDegreeToTest;
    int maximalNumberOfIterations = 300;
    int bestNumberOfIterations=maximalNumberOfIterations;
   // int bestSqareError=INFINITY;
    int valueToReturn = 0;


    vector <float> y; //y ~ ax^2+bx+c
    vector <float> x;

    //check whether there are chromosomes without X and Y:
    vector<ChrCopyNumber>::iterator it;
    int countAutosomes = 0;
    for ( it=chrCopyNumber_.begin() ; it != chrCopyNumber_.end(); it++ ) {
        if (! (it->getChromosome().find("X")!=string::npos || it->getChromosome().find("Y")!=string::npos)) {
            countAutosomes++;
        }
    }
    //fill x and y:
    for ( it=chrCopyNumber_.begin() ; it != chrCopyNumber_.end(); it++ ) {

            if (countAutosomes>0 && !(it->getChromosome().find("X")!=string::npos || it->getChromosome().find("Y")!=string::npos) || countAutosomes==0) {
                // if uniqueMatch, do correction to mappability
                if (uniqueMatch) {
                    //use mappabilityProfile_ and correct
                    for (int i = 0; i< it->getLength(); i++) {
                        if ((it->getValueAt(i)>0)&&(it->getMappabilityProfileAt(i)>minMappabilityPerWindow)) {
                            x.push_back(it->getCGprofileAt(i));
                            y.push_back(it->getValueAt(i)/it->getMappabilityProfileAt(i));
                        }
                    }
                } else {
                    //use a threshold, but correct using notN profile
                    if (it->getMappabilityLength()>0) {
                        for (int i = 0; i< it->getLength(); i++) {
                            if ((it->getValueAt(i)>0)&&(it->getMappabilityProfileAt(i)>minMappabilityPerWindow)) {
                                x.push_back(it->getCGprofileAt(i));
                                y.push_back(it->getValueAt(i)/it->getNotNprofileAt(i));
                            }
                        }
                    } else {
                        for (int i = 0; i< it->getLength(); i++) {
                            if ((it->getValueAt(i)>0)&&(it->getNotNprofileAt(i)>minMappabilityPerWindow)) {
                                x.push_back(it->getCGprofileAt(i));
                                y.push_back(it->getValueAt(i)/it->getNotNprofileAt(i));
                            }
                        }
                    }
                }
            }
    }

    for (int degree=minDegreeToTest;degree<=maxDegreeToTest; degree++) {