-3

次のコードは、プログラムが呼び出されるたびに 2000 万回実行されるため、このコードを可能な限り最適化する方法が必要です。私は高度な C++ プログラマーではないので、この素晴らしいコミュニティの助けを求めています。

int WilcoxinRST::work(GeneSet originalGeneSet, vector<string> randomGenes) {
vector<string> geneIDs;
vector<bool> isOriginal;
vector<int> rank;
vector<double> value;
vector<int> score;
int genesPerSet = originalGeneSet.geneCount();
unsigned int totalGenes, tempScore;
/**
 * Fill the first half of the vectors with original gene set data
 */
totalGenes = genesPerSet * 2;
for (int i = 0; i < genesPerSet; i++) {
    geneIDs.push_back(originalGeneSet.getMemberGeneAt(i));
    isOriginal.push_back(true);
    value.push_back(fabs(expressionLevels.getValue(geneIDs[i], statType)));
}
/**
 * Fill the second half with random data
 */
for (unsigned int i = genesPerSet; i < totalGenes; i++) {
    geneIDs.push_back(randomGenes.at(i - genesPerSet));
    isOriginal.push_back(false);
    value.push_back(fabs(expressionLevels.getValue(geneIDs[i], statType)));
}
totalGenes = geneIDs.size();
/**
 * calculate the scores
 */
if (statType == Fold_Change || statType == T_Statistic
        || statType == Z_Statistic) {
    //      Higher value is a winner
    for (unsigned int i = 0; i < totalGenes; i++) {
        tempScore = 0;
        if (!isOriginal[i]) {
            for (int j = 0; j < genesPerSet; j++) {
                if (value.at(i) > value.at(j)) {
                    tempScore++;
                }
            }

        } else {
            for (unsigned int j = genesPerSet; j < totalGenes; j++) {
                if (value.at(i) > value.at(j)) {
                    tempScore++;
                }
            }

        }

        score.push_back(tempScore);
    }

} else if (statType == FDR_PValue || statType == PValue) {
    // Lower value is a winner
    for (unsigned int i = 0; i < totalGenes; i++) {
        tempScore = 0;
        if (!isOriginal[i]) {
            for (int j = 0; j < genesPerSet; j++) {
                if (value.at(i) < value.at(j)) {
                    tempScore++;
                }
            }

        } else {
            for (unsigned int j = genesPerSet; j < totalGenes; j++) {
                if (value.at(i) < value.at(j)) {
                    tempScore++;
                }
            }

        }

        score.push_back(tempScore);
    }

} else {
    cout << endl << "ERROR. Statistic type not defined." << endl;
}

/**
 * calculate Ua, Ub and U
 */
int U_Original = 0, U_Random = 0, U_Final;
for (int j = 0; j < genesPerSet; j++) {
    U_Original += score[j];
}
for (unsigned int j = genesPerSet; j < totalGenes; j++) {
    U_Random += score[j];
}
U_Final = (U_Original < U_Random) ? U_Original : U_Random;

/**
 * calculate z
 */
double Zn, Zd, Z;
Zn = U_Final - ((genesPerSet * genesPerSet) / 2);
Zd = sqrt(
        (double) (((genesPerSet * genesPerSet
                * (genesPerSet + genesPerSet + 1)))) / 12.0);
Z = Zn / Zd;

/**
 * Return 0/1/2
 * 2: p value < 0.01
 * 1: 0.01 < p value < 0.05
 * 0: p value > 0.05
 */
if (fabs(Z) > 2.303)
    return 2;
else if (fabs(Z) > 1.605)
    return 1;
else
    return 0;
}
4

1 に答える 1

2

コードの複雑さは O(N*N) [ genesPerSet=N] です。しかし、値の順序は関係ないという事実を利用して、O(N•log(N)) で順序付けし、O(N) で「スコア」を計算できます。(潜在的に何千倍も速くなる可能性があります)。

また、合計 N*N の比較があります。次にU_Original + U_Random= N*N であり、U_Random を計算する必要がないことを意味します。また、統計 Zn= Umin-N*N/2;[ Umix=min(U_Original,U_Random)],abs(Zn/Zd) のみが N*N/2 を中心に対称である場合。必要なアルゴリズムは 1 つだけです。

1.- 最初に (const) 参照によって引数を取ることができます:

int WilcoxinRST::work(const GeneSet &originalGeneSet, const vector<string> &randomGenes)

2.- ベクターの遺伝子 ID を入力します。でも使わない？なんで？

3.- 反復できるのは 2 回だけです。

4.- 信号値 (プローブ強度?) を 1 つのベクトルにまとめて保持し、別のベクトルを使用して各項目が何であるかを通知します。単純に 2 つのベクトルに保持します。

5.-スコアベクトルは必要ありません。合計のみが必要です。

6.- なぜ 2000 万回？「統計」の安定性またはBStrapを計算していると思います。おそらく、同じオリジナルのGeneSetを何度も使用しているでしょう。別の質問で、この関数を呼び出すコードを投稿して、毎回値のベクトルを作成して並べ替えることができると思います。

これが最初の新しい O(N•log(N)) コードです。

以下はコードのクリーンアップですが、それでも O(N*N) であり、高速ですが一定の係数だけです。

次に、同じコードですが、元のコードとさらにコメントが混在しています。

これをデバッグして、どうだったか教えてください。

#include<vector>
#include<algorithm>

int WilcoxinRST::work(const GeneSet &originalGeneSet , const vector<string>& randomGenes) 
{
    size_t genesPerSet = originalGeneSet.geneCount();
    std::vector<double> valueOri(genesPerSet), valueRnd(genesPerSet);
    /**
     * Fill the valueOri vector with original gene set data, and valueRnd with random data
     */
    for (size_t i = 0; i < genesPerSet; i++) 
    {
      valueOri[i] = std::fabs(expressionLevels.getValue( originalGeneSet.getMemberGeneAt(i) , statType ));
      valueRnd[i] = std::fabs(expressionLevels.getValue( randomGenes.at(i)                  , statType ));
    }
    std::sort(valueOri.begin(),valueOri.end());
    std::sort(valueRnd.begin(),valueRnd.end());

    /**
     * calculate the scores Ua, Ub and U
     */
    long U_Original ;

    if (statType == Fold_Change || statType == T_Statistic  || statType == Z_Statistic 
        statType == FDR_PValue  || statType == PValue ) 
    {
        //      Higher value is a winner
        size_t j=0;
        for (size_t i = 0; i < genesPerSet /*totalGenes*/; i++)      // i   -  2x
        {   
            while(valueOri[i]  > valueRnd[j]) ++j ;
            U_Original += j;
        }
    } else { cout << endl << "ERROR. Statistic type not defined." << endl;  }

    /**
     * calculate z
     */
    double Zn, Zd, Z;
    Zn = U_Original - ((genesPerSet * genesPerSet) / 2);
    Zd = std::sqrt( (double) (((genesPerSet * genesPerSet* (genesPerSet + genesPerSet + 1)))) / 12.0);
    Z = Zn / Zd;

    /**
     * Return 0/1/2
     * 2: p value < 0.01
     * 1: 0.01 < p value < 0.05
     * 0: p value > 0.05
     */
         if (std::fabs(Z) > 2.303)  return 2;
    else if (std::fabs(Z) > 1.605)  return 1;
    else                            return 0;
}

以下はコードのクリーンアップですが、それでも O(N*N) であり、高速ですが一定の係数だけです。

#include<vector>
using namespace std;
class GeneSet ;
class WilcoxinRST;

int WilcoxinRST::work(const GeneSet &originalGeneSet , const vector<string>& randomGenes) 
{
    size_t genesPerSet = originalGeneSet.geneCount();
    vector<double> valueOri(genesPerSet), valueRnd(genesPerSet);
    /**
     * Fill the valueOri vector with original gene set data, and valueRnd with random data
     */
    for (size_t i = 0; i < genesPerSet; i++) 
    {
        valueOri[i]   =  fabs(expressionLevels.getValue( originalGeneSet.getMemberGeneAt(i) , statType ));
        valueRnd[i]   =  fabs(expressionLevels.getValue( randomGenes.at(i)                  , statType ));
    }
    /**
     * calculate the scores Ua, Ub and U
     */
    long U_Original = 0, U_Random = 0, U_Final;

    if (statType == Fold_Change || statType == T_Statistic  || statType == Z_Statistic) 
    {
        //      Higher value is a winner
        for (size_t i = 0; i < genesPerSet /*totalGenes*/; i++)      // i   -  2x
        {   for (size_t j = 0; j < genesPerSet; j++)   
            {   U_Random  +=  (valueRnd[i]  > valueOri[j]);// i en 2 set=Rnd, j in 1 set=Ori. count how many Ori are less than this Rnd 
                U_Original+=  (valueOri[i]  > valueRnd[j]);// i in 1 set=Ori, j in 2 set=Rnd, count how many Rnd are less than this Ori 
            }
        }
    } else 
    if (statType == FDR_PValue || statType == PValue) 
    {
        // Lower value is a winner
        for (size_t i = 0; i < genesPerSet; i++)   
        {   
            for (size_t j = 0; j < genesPerSet; j++)   
            {   U_Random  +=  (valueRnd[i]  < valueOri[j]);// i en 2 set=Rnd, j in 1 set=Ori. count how many Ori are > than this Rnd 
                U_Original+=  (valueOri[i]  < valueRnd[j]);// i in 1 set=Ori, j in 2 set=Rnd, count how many Rnd are > than this Ori 
            }
        }
    } else { cout << endl << "ERROR. Statistic type not defined." << endl;  }


    U_Final = (U_Original < U_Random) ? U_Original : U_Random;

    /**
     * calculate z
     */
    double Zn, Zd, Z;
    Zn = U_Final - ((genesPerSet * genesPerSet) / 2);
    Zd = sqrt(
            (double) (((genesPerSet * genesPerSet
                    * (genesPerSet + genesPerSet + 1)))) / 12.0);
    Z = Zn / Zd;

    /**
     * Return 0/1/2
     * 2: p value < 0.01
     * 1: 0.01 < p value < 0.05
     * 0: p value > 0.05
     */
         if (fabs(Z) > 2.303)       return 2;
    else if (fabs(Z) > 1.605)       return 1;
    else                            return 0;
}

同じコードですが、元のコードとより多くのコメントが混在しています。

int WilcoxinRST::work(const GeneSet &originalGeneSet , const vector<string>& randomGenes) 
{
    size_t genesPerSet = originalGeneSet.geneCount();
    unsigned int totalGenes, tempScore;
    totalGenes = genesPerSet * 2;

    //vector<string> geneIDs;
    //vector<bool>   isOriginal;
    //vector<int>    rank;
    vector<double> valueOri(genesPerSet), valueRnd(genesPerSet);
    //vector<int>    score;
    /**
     * Fill the first half of the vectors with original gene set data
     */

    for (size_t i = 0; i < genesPerSet; i++) 
    {
        //geneIDs.push_back( originalGeneSet.getMemberGeneAt(i)  );
        //isOriginal.push_back(true);

        valueOri[i]   =  fabs(expressionLevels.getValue( originalGeneSet.getMemberGeneAt(i) , statType ));
        valueRnd[i]   =  fabs(expressionLevels.getValue( randomGenes.at(i)                  , statType ));
    }
    /**
     * Fill the second half with random data
     */
    //for (unsigned int i = genesPerSet; i < totalGenes; i++) {
    //  geneIDs.push_back(randomGenes.at(i - genesPerSet));
    //  isOriginal.push_back(false);
    //  value.push_back(fabs(expressionLevels.getValue(geneIDs[i], statType)));
    //}
    //totalGenes = geneIDs.size();
    /**
     * calculate the scores
     */
        /**
     * calculate Ua, Ub and U
     */
    long U_Original = 0, U_Random = 0, U_Final;
    //for (int j = 0; j < genesPerSet; j++) // j in 1 set=Ori. count how many Ori are less than this Rnd
    //{
    //  U_Original += score[j];
    //}
    //for (unsigned int j = genesPerSet; j < totalGenes; j++) // j in 2 set=Rnd, count how many Rnd are less than this Ori 
    //{
    //  U_Random += score[j];
    //}

    if (statType == Fold_Change || statType == T_Statistic  || statType == Z_Statistic) 
    {
        //      Higher value is a winner
        for (size_t i = 0; i < genesPerSet /*totalGenes*/; i++)      // i   -  2x
        {   //tempScore = 0;
            //if (!isOriginal[i])  // i en 2 set=Rnd, j in 1 set=Ori. count how many Ori are less than this Rnd 
            for (size_t j = 0; j < genesPerSet; j++)   
            {   U_Random  +=  (valueRnd[i]  > valueOri[j]);// i en 2 set=Rnd, j in 1 set=Ori. count how many Ori are less than this Rnd 
                U_Original+=  (valueOri[i]  > valueRnd[j]);// i in 1 set=Ori, j in 2 set=Rnd, count how many Rnd are less than this Ori 
            }
            //} else 
            //{
            //  for (unsigned int j = genesPerSet; j < totalGenes; j++)  // i in 1 set=Ori, j in 2 set=Rnd, count how many Rnd are less than this Ori 
            //  {   if (value.at(i) > value.at(j)) {    tempScore++;        }
            //  }

            //}
            //score.push_back(tempScore);
        }

    } else 
    if (statType == FDR_PValue || statType == PValue) 
    {
        // Lower value is a winner
        for (size_t i = 0; i < genesPerSet; i++)   
        {   
            for (size_t j = 0; j < genesPerSet; j++)   
            {   U_Random  +=  (valueRnd[i]  < valueOri[j]);// i en 2 set=Rnd, j in 1 set=Ori. count how many Ori are > than this Rnd 
                U_Original+=  (valueOri[i]  < valueRnd[j]);// i in 1 set=Ori, j in 2 set=Rnd, count how many Rnd are > than this Ori 
            }
            //} else 
            //{
            //  for (unsigned int j = genesPerSet; j < totalGenes; j++)  // i in 1 set=Ori, j in 2 set=Rnd, count how many Rnd are less than this Ori 
            //  {   if (value.at(i) > value.at(j)) {    tempScore++;        }
            //  }

            //}
            //score.push_back(tempScore);
        }

        //for (unsigned int i = 0; i < totalGenes; i++) 
        //{   tempScore = 0;
        //  if (!isOriginal[i]) 
        //  {   for (int j = 0; j < genesPerSet; j++) {
        //          if (value.at(i) < value.at(j)) {  // Rnd i < Ori j increm U_Random
        //              tempScore++;
        //          }
        //      }

        //  } else {
        //      for (unsigned int j = genesPerSet; j < totalGenes; j++) { // Ori i < Rnd j. Increm U_Original
        //          if (value.at(i) < value.at(j)) {
        //              tempScore++;
        //          }
        //      }

        //  }

        //  score.push_back(tempScore);
        //}

    } else { cout << endl << "ERROR. Statistic type not defined." << endl;  }


    U_Final = (U_Original < U_Random) ? U_Original : U_Random;

    /**
     * calculate z
     */
    double Zn, Zd, Z;
    Zn = U_Final - ((genesPerSet * genesPerSet) / 2);
    Zd = sqrt(
            (double) (((genesPerSet * genesPerSet
                    * (genesPerSet + genesPerSet + 1)))) / 12.0);
    Z = Zn / Zd;

    /**
     * Return 0/1/2
     * 2: p value < 0.01
     * 1: 0.01 < p value < 0.05
     * 0: p value > 0.05
     */
    if (fabs(Z) > 2.303)
        return 2;
    else if (fabs(Z) > 1.605)
        return 1;
    else
        return 0;
}
于 2013-03-28T07:55:26.023 に答える