Skip to content

Unmapped.out.mate* names and sequences as STAR writes them - #364

Open
BenjaminDEMAILLE wants to merge 1 commit into
fix/input-output-smallfrom
fix/unmapped-mate-names
Open

BenjaminDEMAILLE wants to merge 1 commit into
fix/input-output-smallfrom
fix/unmapped-mate-names

Conversation

@BenjaminDEMAILLE

Copy link
Copy Markdown
Contributor

Part of #354 (the Unmapped.out.mate* name difference reported in #358; it has no finding number in the audit report).

Stacked on #358, because this touches the same FASTQ code and uses the FASTA reads it adds. Only the last commit is this PR. Its base branch here is fix/input-output-small, and #358 is itself stacked on #353 and #352.

Problem

--outReadsUnmapped Fastx wrote @q9 where STAR writes @q9 0:N: , so the files differed on every record. Measuring this turned up three more differences in the same files:

  • A FASTQ mate 2 record used mate 2's own name, where STAR uses mate 1's.
  • Half-mapped pairs were written clipped. STAR writes the mates as read, and the other unmapped path in the same function already did.
  • The filter flag and mate status were not written at all.

What STAR does

Observed output of Homebrew STAR 2.7.11b, source not read.

single-end   <m><name> 0:<flag>: 
pair, mate 1 <m><name> 0:<flag>:  <mate1 mapped><mate2 mapped>
pair, mate 2 <m><name> 1:<flag>:  <mate1 mapped><mate2 mapped>

<m> is @ or >. The text between the flag's : and the status is empty in every run I could make, so there are two spaces before the status in pairs, and a trailing space in single-end. The audit calls that the barcode field. I could not confirm that, see "Not matched".

Part Observed
Name Cut at --readNameSeparator as elsewhere. With --outSAMreadID Number, the number, on both mates
Mate 2's name, FASTQ Mate 1's: @d.0.1 1:N: 00 for d.0.1 / d.0.2
Mate 2's name, FASTA Its own: >d.0.2 1:N: 00
Flag, FASTQ Y if the first word after the name has :Y: at its second byte, else N
Flag, FASTA Always N, even for >c1 1:Y:0
Flag, pairs Mate 1's on both records. 1:Y:0 / 2:N:0 gives Y on both, and 1:N:0 / 2:Y:0 gives N on both
Status 10 mate 1 mapped, 01 mate 2 mapped, 00 otherwise, including a pair over --outFilterMultimapNmax
Sequence and quality As read, not clipped. The + line is bare, also for +q14

Flag on single-end FASTQ, as written by STAR:

Name line Flag
@a1 1:Y:0:ACGT, @a3 2:Y:18:GG, @a12 1:Y:0 2:N:0, @a13 1:Y:, @a9/1/x y:Y: Y
@a8<TAB>b:Y:c Y
@a2 1:N:0:ACGT, @a4 Y, @a5 :Y:, @a6 1:Y, @a10 ab:Y:, @a11 1:y:0 N
@a7 x<TAB>y:Y:z N

Change

  • EncodedRead.filter_flag, set by the new filter_flag for FASTQ and N for FASTA. Reading it costs one short scan per read and a byte in the struct. Its only use is here.
  • unmapped_name builds the name, and the four call sites in src/lib.rs use it. UnmappedFastqWriter itself did not need to change.
  • Half-mapped pairs write paired_read.mate* rather than the clipped copies.

Other EncodedRead constructors are test code in src/solo/mod.rs; each gains one line.

Measured against STAR

cmp on Unmapped.out.mate*. Exit codes were 0 on both tools unless shown.

Run Records per file Result
SE FASTQ, comments, /1, tab, +q14 8 identical
SE FASTA 3 identical
SE FASTQ, flag spellings above 12 + 2 identical
SE FASTA with 1:Y:0 2 identical
Pair, FASTQ, flag on one mate each 2 identical, both mates
Pair, FASTQ, 10 name shapes 9 identical, both mates
Pair, FASTA, same names 9 identical, both mates
--readNameSeparator _ (pair); q and - (SE) 9; 8; 8 identical
--outSAMreadID Number, SE and pair 8; 9 identical
Half-mapped pair, 10, with --clip5pNbases 3 3 --clip3pNbases 2 2 2 identical, both mates
Half-mapped pair, 01, same clipping 2 identical, both mates
Half-mapped, 10 and 01, FASTA 2 identical, both mates
Half-mapped, --outSAMreadID Number 2 identical, both mates
Pair over --outFilterMultimapNmax 1, both mates mapped 2 identical, both mates
Strict filters with clipping, pair 4 identical, both mates
--outFilterType BySJout, pair 9 identical, both mates
Half-mapped pair over --outFilterMultimapNmax 1 3 differs, one name, see below
Mixed FASTQ and FASTA pair 0 both empty; exit 104 and 1

I did not re-run #358's own binary for a before column. Its description records >junk against STAR's >junk 0:N: , and main gave @q9 against @q9 0:N: here.

Not matched

  • A half-mapped pair over --outFilterMultimapNmax: STAR writes @multi_polyA 0:N: 00, this writes 10. On this branch a half-mapped pair always has one result, so it is never seen as over the limit. That is finding 2 of Tracking: fixes from the whole-repo audit (main at 32fbbc1) #354, and I left it there.
  • Barcode field: never observed non-empty. I expect STARsolo to use it, but on macOS STAR stops on --soloType runs with could not open input file /geneInfo.tab. So I do not know what STAR writes there.
  • STARsolo runs: the solo read loops here write nothing to these files, on main too. Not changed and not compared.
  • Mixed pairs (FASTQ with FASTA): both tools refuse them with no records, but STAR exits 104 and this exits 1.

Deliberate divergences

None added.

Tests

Each failed with the fix taken out, and passes with it:

  • tests/alignment_features.rs: test_unmapped_mate_files_match_star. Compares the whole of each file: single-end FASTQ and FASTA, a FASTQ and a FASTA pair, and half-mapped pairs with clipping, with both 10 and 01. I ran STAR on the same genome and reads and its output is byte-identical to the test's expected text.
  • src/io/fastq.rs: filter_flag_needs_y_as_the_second_colon_field_of_a_one_byte_first_field.

cargo test --release, cargo clippy --all-targets and cargo fmt --check are clean on macOS.

Not checked

  • test/parity_suite.sh was not run. This changes those files only, so the SAM output comparison harnesses were not run either.
  • Any tool that reads these files, and depends on the old names, was not checked. With @q9 0:N: , a tool taking the whole line as the name would see the extra fields, as it does with STAR.
  • Only macOS was run.

🤖 Generated with Claude Code

`@q9` becomes `@q9 0:N: `, and pairs get `@p 0:N:  10` / `@p 1:N:  10`: the
0-based mate, the FASTQ filter flag, an empty field, and whether each mate
mapped. A FASTQ mate 2 record carries mate 1's name and flag; a FASTA one keeps
its own name. Half-mapped pairs are written as read, not clipped.

Format observed in STAR 2.7.11b's output, source not read.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>

This branch has not been deployed

No deployments
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant