getzlab/rnaseqc

5'/3' bias edge case bugs

Closed this issue · 1 comments

A few bugs in calculation of 5'/3' bias, specifically for the edge case where only 1 gene is included in the calculation.

  1. In calculation of 25th and 75th percentiles, the index to read at is not checked against the size of the ratios vector. If the ratios vector is length 1, this leads to reading from past the end of the vector, giving undefined behavior. Also, the behaviors in both the if and else blocks are currently mathematically identical, I believe.

rnaseqc/src/RNASeQC.cpp

Lines 500 to 521 in 19967f3

double index = .25 * ratios.size();
if (index > floor(index))
{
index = ceil(index);
ratio25 = ratios[static_cast<int>(index)];
}
else
{
index = ceil(index);
ratio25 = (ratios[static_cast<int>(index)] + ratios[static_cast<int>(index)])/2.0;
}
index = .75 * ratios.size();
if (index > floor(index))
{
index = ceil(index);
ratio75 = ratios[static_cast<int>(index)];
}
else
{
index = ceil(index);
ratio75 = (ratios[static_cast<int>(index)] + ratios[static_cast<int>(index)])/2.0;
}

  1. In the code to compute medians, the edge case when size is 1 hits a similar bug.

rnaseqc/src/Metrics.h

Lines 152 to 156 in 19967f3

if (size % 2)
{
double value = static_cast<double>(*(iterator++));
return (value + static_cast<double>(*iterator)) / 2.0;
}

Thanks for pointing these out. Good catches!