PSI=0 when 0 reads

[mirimia]

Hi guys,

I haven't used this for a while, but now that I'm back to it I think I saw a relatively important problem (although probably didn't really have any effect on anything yet): when there are 0 reads (for inc and excl) it seems the PSI is 0. This can be a big problem if no minimum coverage is taken into account, as sometime we've done for instance for the IR paper. Is that correct? In my tables it used to be an NA when that happened.

Thanks, Manu

[mirimia]

Or even worse, these seem to have PSI=100 with 0 reads?!? What am I missing?!

100 N,N,N,Bn,S@0,0 0 N,N,N,Bn,S@0,0 100 N,N,N,Bn,S@0,0

[kcha]

Hi Manuel,

Can you give us a complete example that we can perhaps try to replicate on our end?

AFAIK, aside from Tim's addition of the counts to the quality score, we haven't changed anything in the PSI formulas or its output. I could be wrong...

Kevin

[timbitz]

Yah, I agree with Kevin, the psi score should be as-is, though that doesn't mean I didn't accidentally fiddle with it somehow..

I can try and look into it, though I have to admit I didn't always have the best picture as to what was going on there to begin with... ;-)

On Oct 21, 2014, at 10:21 PM, Kevin Ha notifications@github.com<mailto:notifications@github.com> wrote:

Hi Manuel,

Can you give us a complete example that we can perhaps try to replicate on our end?

AFAIK, aside from Tim's addition of the counts to the quality score, we haven't changed anything in the PSI formulas or its output. I could be wrong...

Kevin

— Reply to this email directly or view it on GitHubhttps://github.com/vastgroup/vast-tools/issues/29#issuecomment-60027909.

[mirimia]

That's weird... OK, I'll have a look at it, as maybe it missed something. \@0,0 is directly from the read count? I see when there are reads in either inc or exc, but not the other, the "0" count is 0.00, right?

Thanks, Manu

[mirimia]

Any chance when there is only 1 read the \@inc,exc is \@inc,exc?

Thanks,

[mirimia]

OK, I think I got it. It requires $totalN < 2. Sorry about it.

my $totalN = $total_reads; # no point in re-assigning really. my($pPSI, $exValOfInc, $exValOfExc) = (0, 0, 0); unless($PSI_complex eq "NA" or $totalN < 2) { $pPSI = $PSI_complex / 100;

$exValOfInc = $pPSI * $totalN;

$exValOfExc = (1-$pPSI) * $totalN;

$exValOfInc = sprintf("%.2f", $pPSI * $totalN); $exValOfExc = sprintf("%.2f", (1-$pPSI) * $totalN); }

ALTER QUAL OUTPUT HERE>>

$Q .= "\@$exValOfInc,$exValOfExc";