Skip to content
Merged
Show file tree
Hide file tree
Changes from 3 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
112 changes: 111 additions & 1 deletion pwiz_tools/Osprey/Osprey.Tasks/Pass2FdrSidecar.cs
Original file line number Diff line number Diff line change
Expand Up @@ -195,6 +195,24 @@ internal static FeatureContributions ComputeAndPersist(
"them here from the reconciled features (reused distributed-run code path).",
missingPass2, totalFiles));
var swPass2 = Stopwatch.StartNew();

// Stage 6's post-rescore overlay calls FdrEntry.ResetScores(), which clears
// eight fields. Pass 2 recomputes only five of them, and which five depends
// on the mode and on whether a survivor is on-stratum, so three can reach the
// 2nd-pass sidecar at their reset defaults (issue #4553):
// Score - no frozen mode wrote one back at all
// Pep - written only on the on-stratum path
// RunProteinQvalue - written by NO mode; first-pass protein FDR is its
// only producer, and the second-pass one writes
// ExperimentProteinQvalue instead
// Seeding all three from the 1st-pass sidecar reproduces exactly what the
// distributed route has in hand at this point, which is why that route never
// showed the loss: it must rehydrate from that same sidecar, and the sidecar
// carries all seven scalars. Whatever pass 2 genuinely recomputes then
// overwrites the seed. Done ahead of the mode dispatch because the loss is
// not specific to one mode.
RestorePass1Scalars(ctx, perFileEntries, inputByFileName);

// --model-diagnostics needs the resident 2nd-pass model: its feature
// contributions feed the pass-2 model view, and the projection 2nd pass
// streams through a sink and produces none. Route --model-diagnostics to
Expand Down Expand Up @@ -506,6 +524,92 @@ void FlushPass2File(string fileName, IReadOnlyList<FdrScoreRecord> records)
return pass2Contributions;
}

/// <summary>
/// Re-seed each survivor's <see cref="FdrEntry.Score"/>, <see cref="FdrEntry.Pep"/> and
/// <see cref="FdrEntry.RunProteinQvalue"/> from that file's
/// <c>.1st-pass.fdr_scores.bin</c>.
///
/// <para>These are the three of <see cref="FdrEntry.ResetScores"/>'s eight fields that
/// pass 2 does not reliably recompute: no frozen mode wrote <c>Score</c> at all,
/// <c>Pep</c> is written only for on-stratum survivors, and <c>RunProteinQvalue</c> is
/// written by no mode at all. Left unseeded they reach the 2nd-pass sidecar at their
/// reset defaults, where a q-value of 1.0 reads as a confident rejection and a
/// <c>Score</c> of 0 sits exactly ON the discriminant's accept/reject boundary
/// (issue #4553).</para>
///
/// <para>Seeding, not overriding: whatever pass 2 genuinely recomputes is written after
/// this and wins. What is left is the pass-1 value, which is precisely what the
/// distributed route holds at the same point - it rehydrates from this same sidecar -
/// so the two routes agree by construction rather than by coincidence.</para>
///
/// <para>An entry Stage 6 did not touch already holds these values, so the write is a
/// no-op for it; a gap-fill entry is absent from the sidecar and correctly keeps the
/// defaults, which is where the distributed route leaves it too (its own overlay runs
/// before gap-fill appends). One file's records stream at a time.</para>
/// </summary>
private static void RestorePass1Scalars(
PipelineContext ctx,
List<KeyValuePair<string, List<FdrEntry>>> perFileEntries,
IReadOnlyDictionary<string, string> inputByFileName)
{
int nRestored = 0;
int filesRead = 0;
var unreadable = new List<string>();
foreach (var kvp in perFileEntries)
{
if (!inputByFileName.TryGetValue(kvp.Key, out string inputFile))
continue;
string pass1Path = FdrScoresSidecar.Pass1Path(inputFile);
if (!File.Exists(pass1Path))
{
unreadable.Add(kvp.Key);
continue;
}
var byEntryId = new Dictionary<uint, FdrEntry>(kvp.Value.Count);
foreach (var e in kvp.Value)
byEntryId[e.EntryId] = e;

int nFile = 0;
bool ok = FdrScoresSidecar.ReadRecords(
pass1Path, FdrScoresSidecar.Pass.FirstPass,
rec =>
{
if (!byEntryId.TryGetValue(rec.EntryId, out FdrEntry entry))
return;
entry.Score = rec.Score;
entry.Pep = rec.Pep;
entry.RunProteinQvalue = rec.RunProteinQvalue;
nFile++;
});
if (ok)
{
filesRead++;
nRestored += nFile;
}
else
{
unreadable.Add(kvp.Key);
}
}

// A missing or corrupt 1st-pass sidecar is reported, not thrown on: none of these
// three fields is an input to a pass-2 computation, so the run stays correct in
// every other respect, and the frozen modes have their own fail-fast for a sidecar
// they genuinely need. Silence, though, would leave the reset defaults looking
// computed.
if (unreadable.Count > 0)
{
ctx.LogWarning(string.Format(
"1st-pass Score/Pep/RunProteinQvalue could not be restored for {0} file(s) " +
"(no readable 1st-pass sidecar): [{1}]. Their 2nd-pass sidecars will carry " +
"reset defaults for peaks Stage 6 changed.",
unreadable.Count, string.Join(", ", unreadable)));
}
ctx.LogVerbose(string.Format(
"Restored 1st-pass Score/Pep/RunProteinQvalue onto {0} survivor(s) across {1} file(s).",
nRestored, filesRead));
}

/// <summary>
/// OSPREY_PASS2_QVALUE=transfer-compete (full-population form). Recompute the reported
/// precursor q-values + PEP by re-running the target-decoy competition over the ENTIRE
Expand Down Expand Up @@ -769,7 +873,13 @@ private static bool ComputePass2TransferCompeteFull(
(e.EntryId, e.Charge, e.ScanNumber), out double[] feats) &&
feats != null && feats.Length == nFeatures)
{
fileScores[e.EntryId] = scorer.Score(feats);
double frozenScore = scorer.Score(feats);
fileScores[e.EntryId] = frozenScore;
// This is the score the entry COMPETES on below, so it is the one
// the 2nd-pass sidecar must carry. RestorePass1Scalars seeded the
// 1st-pass value, which is what a survivor whose features did not
// resolve keeps - and which is what it competes on too.
e.Score = frozenScore;
}
}
// featByIdentity released here (one file resident at a time).
Expand Down
234 changes: 234 additions & 0 deletions pwiz_tools/Osprey/Regression/FdrSidecars.ps1
Original file line number Diff line number Diff line change
@@ -0,0 +1,234 @@
<#
Per-file FDR score sidecar comparison.

Two consumers, one decoder:
* regression.ps1's four-task-chain leg -- distributed route vs straight-through, C# only.
* ai/scripts/Osprey/Compare/Compare-FdrSidecars-Crossimpl.ps1 -- C# vs Rust.

The four-task chain leg asserted only Compare-BlibFull, and the blib carries no protein
q-value, so a route that writes a different RunProteinQvalue into every
<stem>.2nd-pass.fdr_scores.bin passed green (issue #4553: 32,450 of 260,419 records differ
on StellarGenDecoyEntrap, 1.57% at 82 files). Peptide counts, protein-group counts and the
blib are all identical while it happens, so nothing the gate already reads can see it.

The cross-impl gate had the same blind spot from the other direction: it compares the
Stage 7 protein FDR dump (per-protein-GROUP columns) and the blib, neither of which carries
a per-entry SVM score or run protein q. Both implementations dropped the same two fields,
so they agreed on the wrong value and nothing was red.

This compares the sidecars themselves, which is where the distributed tasks' per-file
output actually lands.

Record layout (Osprey.IO\FdrScoresSidecar.cs): 32-byte header, 60-byte records,
entry_id u32 @0, score f64 @4, run_precursor_q @12, run_peptide_q @20,
experiment_precursor_q @28, experiment_peptide_q @36, pep @44, run_protein_q @52.

The decode + compare runs as compiled C#, not PowerShell. A per-record PowerShell loop
took over 10 minutes on one Astral 3-file pass (6.2M records) and would be unusable at the
82-file scale these gates are meant to reach; the compiled form is seconds.
#>

if (-not ([System.Management.Automation.PSTypeName]'OspreyFdrSidecarComparer').Type) {
Add-Type -TypeDefinition @'
using System;
using System.Collections.Generic;
using System.IO;

public class FdrSidecarDiff
{
public bool ReadableExpected;
public bool ReadableActual;
public string ProblemExpected;
public string ProblemActual;
public long CountExpected;
public long CountActual;
public long Compared;
public long[] FieldCounts;
public string[] FirstExample;
}

public static class OspreyFdrSidecarComparer
{
// Byte offsets within the 60-byte record, in report order.
private static readonly int[] Offsets = { 4, 12, 20, 28, 36, 44, 52 };
private const int HeaderLen = 32;
private const int RecordLen = 60;
private const byte ExpectedVersion = 3;
private static readonly byte[] Magic = { 0x4F, 0x53, 0x50, 0x52, 0x59, 0x46, 0x44, 0x52 }; // OSPRYFDR

public static string[] FieldNames = {
"score", "run_precursor_qvalue", "run_peptide_qvalue",
"experiment_precursor_qvalue", "experiment_peptide_qvalue",
"pep", "run_protein_qvalue"
};

public static FdrSidecarDiff Compare(string pathExpected, string pathActual, double tolerance)
{
var result = new FdrSidecarDiff
{
FieldCounts = new long[Offsets.Length],
FirstExample = new string[Offsets.Length],
};

byte[] a = ReadIfValid(pathExpected, out long na, out string problemA);
byte[] b = ReadIfValid(pathActual, out long nb, out string problemB);
result.ReadableExpected = a != null;
result.ReadableActual = b != null;
result.ProblemExpected = problemA;
result.ProblemActual = problemB;
result.CountExpected = na;
result.CountActual = nb;
if (a == null || b == null)
return result;

// entry_id -> record offset for the actual side, then walk the expected side.
var offsetById = new Dictionary<uint, int>((int)nb);
for (long i = 0; i < nb; i++)
{
int off = HeaderLen + (int)(i * RecordLen);
offsetById[BitConverter.ToUInt32(b, off)] = off;
}

for (long i = 0; i < na; i++)
{
int offA = HeaderLen + (int)(i * RecordLen);
uint entryId = BitConverter.ToUInt32(a, offA);
int offB;
if (!offsetById.TryGetValue(entryId, out offB))
continue;
result.Compared++;
for (int f = 0; f < Offsets.Length; f++)
{
double va = BitConverter.ToDouble(a, offA + Offsets[f]);
double vb = BitConverter.ToDouble(b, offB + Offsets[f]);
if (Math.Abs(va - vb) <= tolerance)
continue;
result.FieldCounts[f]++;
if (result.FirstExample[f] == null)
{
result.FirstExample[f] = string.Format(
"entry_id={0} {1:R} -> {2:R}", entryId, va, vb);
}
}
}
return result;
}

/// Read a sidecar whose header validates, else null with a reason.
///
/// Magic and version are checked, not just the size invariant. A newer writer that grows
/// the record (the v4 68-byte layout on the #4522 branch adds experiment_aggregate_score
/// at [60..68]) must be REFUSED by name here rather than decoded at the old stride: the
/// size check alone rejects it only because 68 and 60 happen not to divide alike, which
/// is luck, not a guard. Silently misreading a field offset is the failure mode this
/// comparison exists to catch, so it must not be one this comparison can commit.
private static byte[] ReadIfValid(string path, out long count, out string problem)
{
count = 0;
problem = null;
if (!File.Exists(path))
{
problem = "file not found";
return null;
}
byte[] data = File.ReadAllBytes(path);
if (data.Length >= 9)
{
for (int i = 0; i < Magic.Length; i++)
{
if (data[i] != Magic[i])
{
problem = "wrong magic bytes (not an Osprey FDR sidecar)";
return null;
}
}
if (data[8] != ExpectedVersion)
{
problem = string.Format(
"sidecar format version {0}, but this comparison decodes version {1} " +
"({2}-byte records). Update FdrSidecars.ps1 for the newer layout.",
data[8], ExpectedVersion, RecordLen);
return null;
}
}
if (data.Length < HeaderLen)
return null;
ulong n = BitConverter.ToUInt64(data, 16);
if ((ulong)data.Length != (ulong)HeaderLen + n * RecordLen)
return null;
count = (long)n;
return data;
Comment thread
brendanx67 marked this conversation as resolved.
Outdated
}
}
'@
}

function Compare-FdrSidecars {
<#
Compare every <stem>.<Pass>-pass.fdr_scores.bin between two run directories, all seven
scalar fields at Tolerance, matched by stem and then by entry_id. Returns Pass + Issues.

-Pass selects which sidecar: 2 (default) is the one #4553 is about; 1 is the
pre-reconciliation write, compared by the cross-impl gate so a 2nd-pass failure can be
read as "pass 2 dropped it" rather than "the two runs diverged upstream".
#>
param(
[Parameter(Mandatory = $true)][string]$ExpectedDir,
[Parameter(Mandatory = $true)][string]$ActualDir,
[ValidateSet(1, 2)][int]$Pass = 2,
[double]$Tolerance = 1e-9
)

$issues = [System.Collections.Generic.List[string]]::new()
$suffix = if ($Pass -eq 1) { '.1st-pass.fdr_scores.bin' } else { '.2nd-pass.fdr_scores.bin' }
$expected = @{}
foreach ($f in Get-ChildItem -File -Path $ExpectedDir -Filter "*$suffix" -ErrorAction SilentlyContinue) {
$expected[$f.Name.Substring(0, $f.Name.Length - $suffix.Length)] = $f.FullName
}
$actual = @{}
foreach ($f in Get-ChildItem -File -Path $ActualDir -Filter "*$suffix" -ErrorAction SilentlyContinue) {
$actual[$f.Name.Substring(0, $f.Name.Length - $suffix.Length)] = $f.FullName
}

if ($expected.Count -eq 0) {
$issues.Add("no $suffix files in expected dir $ExpectedDir")
return @{ Pass = $false; Issues = $issues; Compared = 0 }
}
foreach ($stem in $expected.Keys) {
if (-not $actual.ContainsKey($stem)) { $issues.Add("missing $stem$suffix in $ActualDir") }
}
foreach ($stem in $actual.Keys) {
if (-not $expected.ContainsKey($stem)) { $issues.Add("unexpected $stem$suffix in $ActualDir") }
}

$nCompared = 0
foreach ($stem in ($expected.Keys | Sort-Object)) {
if (-not $actual.ContainsKey($stem)) { continue }

$diff = [OspreyFdrSidecarComparer]::Compare($expected[$stem], $actual[$stem], $Tolerance)
if (-not $diff.ReadableExpected -or -not $diff.ReadableActual) {
# Name the reason. "Not readable" on a version bump reads as corruption and sends
# the next person looking at the wrong thing.
$why = if (-not $diff.ReadableExpected) { $diff.ProblemExpected } else { $diff.ProblemActual }
$side = if (-not $diff.ReadableExpected) { 'expected' } else { 'actual' }
$issues.Add("$stem$suffix unreadable on the $side side: $why")
continue
}
if ($diff.CountExpected -ne $diff.CountActual) {
$issues.Add("$stem$suffix record count $($diff.CountExpected) -> $($diff.CountActual)")
}
$nCompared += $diff.Compared

Comment thread
brendanx67 marked this conversation as resolved.
# Per-field tallies so the summary names WHICH field drifted, not just that
# something did -- the fields have very different failure meanings.
for ($f = 0; $f -lt [OspreyFdrSidecarComparer]::FieldNames.Length; $f++) {
if ($diff.FieldCounts[$f] -gt 0) {
$issues.Add(("{0}: {1} differs on {2} record(s); first {3}" -f
$stem, [OspreyFdrSidecarComparer]::FieldNames[$f],
$diff.FieldCounts[$f], $diff.FirstExample[$f]))
}
}
}

return @{ Pass = ($issues.Count -eq 0); Issues = $issues; Compared = $nCompared }
}
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@ RetentionTimes <rows> 353298
RetentionTimes retentionTime 3446613.4058412956 1.961938983333 23.813372533333
RetentionTimes startTime 3914888.547512469 1.202014933333 23.7654244
RetentionTimes endTime 3962646.571064724 2.0269066 23.870514
RetentionTimes score 8287.14162545012 3.0744635061181826E-05 1
RetentionTimes score 8287.141625450126 3.0744635061181826E-05 1
OspreyRunScores <rows> 117783
OspreyRunScores RunQValue 66.91619461292912 3.0744635061181826E-05 0.009958467928011075
OspreyRunScores DiscriminantScore 0 0 0
Expand Down
Loading