Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Test secondary alignments #365

Merged
merged 2 commits into from
Nov 17, 2023
Merged

Test secondary alignments #365

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
fd798bd
Use an asterisk for empty SEQ fields in SAM
marcelm Nov 17, 2023
35bbd64
Test secondary alignments
marcelm Nov 17, 2023
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
10 changes: 5 additions & 5 deletions src/sam.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -83,7 +83,7 @@ void Sam::add_unmapped(const KSeq& record, uint16_t flags) {
sam_string.append("\t");
sam_string.append(std::to_string(flags));
sam_string.append("\t*\t0\t" SAM_UNMAPPED_MAPQ_STRING "\t*\t*\t0\t0\t");
sam_string.append(record.seq);
append_seq(record.seq);
append_qual(record.qual);
append_rg_and_newline();
}
Expand All @@ -105,7 +105,7 @@ void Sam::add_unmapped_mate(const KSeq& record, uint16_t flags, const std::strin
sam_string.append("\t");
sam_string.append(std::to_string(mate_pos + 1));
sam_string.append("\t0\t");
sam_string.append(record.seq);
append_seq(record.seq);
append_qual(record.qual);
append_rg_and_newline();
}
Expand Down Expand Up @@ -176,11 +176,11 @@ void Sam::add_record(
sam_string.append("\t");

if (flags & SECONDARY) {
sam_string.append("*");
append_seq("");
} else if (flags & REVERSE) {
sam_string.append(query_sequence_rc);
append_seq(query_sequence_rc);
} else {
sam_string.append(query_sequence);
append_seq(query_sequence);
}

if (!(flags & UNMAP)) {
Expand Down
5 changes: 5 additions & 0 deletions src/sam.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -86,6 +86,11 @@ class Sam {

private:
void add_record(const std::string& query_name, uint16_t flags, const std::string& reference_name, uint32_t pos, uint8_t mapq, const Cigar& cigar, const std::string& mate_reference_name, uint32_t mate_pos, int32_t template_len, const std::string& query_sequence, const std::string& query_sequence_rc, const std::string& qual, int ed, int aln_score, const Details& details);

void append_seq(const std::string& seq) {
sam_string.append(seq.empty() ? "*" : seq);
}

void append_qual(const std::string& qual) {
sam_string.append("\t");
sam_string.append(qual.empty() ? "*" : qual);
Expand Down
4 changes: 2 additions & 2 deletions tests/phix.pe.sam
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -89,5 +89,5 @@ rescuable.43 83 NC_001422.1 3137 60 120S14=1X4=1X9=1X4=1X4=1X9=1X4=1X4=1X4=1X9=1
rescuable.43 163 NC_001422.1 2955 60 1X177=1X122= = 3137 362 NTTGCTTCTGCTCTTGCTGGTGGCGCCATGTCTAAATTGTTTGGAGGCGGTCAAAAAGCCGCCTCCGGTGGCATTCAAGGTGATGTGCTTGCTACCGATAACAATACTGTAGGCATGGGTGATGCTGGTATTAAATCTGCCATTCAAGGCTCTAATGTTCCTAACCCTGATGAGGCCGTCCCTAGTTTTGTTTCTGGTGCTATGGCTAAAGCTGGTAAAGGACTTCTTGAAGGTACGTTGCAGGCTGGCACTTCTGCCGTTTCTGATAAGTTGCTTGATTTGGTTGGACTTGGTGGCAAGT #8BCCGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGFEFGGGCFGGGGFGGGGGGGGGGGGGGGGFGGGFFCGGGFGGGGGGGGGGGGGGGFGGDEG>C9:FCFFDFGCFGGGGGGFFGFFFFF4AAFFFFFFFFFFFFFBEFFFFFE=@@8C@3>BF=4=>;?6?@F@;E<CE@FF<E3.<?37);?E?94<>:<A>BFF?F6?F?)))<?FFF### NM:i:2 AS:i:602 RG:Z:1
not.rescuable 89 NC_001422.1 3017 60 116=1X183=1X = 3017 0 CTCCGGTGGCATTCAAGGTGATGTGCTTGCTACCGATAACAATACTGTAGGCATGGGTGATGCTGGTATTAAATCTGCCATTCAAGGCTCTAATGTTCCTAACCCTGATGAGGCCGTCCCTAGTTTTGTTTCTGGTGCTATGGCTAAAGCTGGTAAAGGACTTCTTGAAGGTACGTTGCAGGCTGGCACTTCTGCCGTTTCTGATAAGTTGCTTGATTTGGTTGGACTTGGTGGCAAGTCTGCCGCTGATAAAGGAAAGGATACTCGTGATTATCTTGCTGCTGCATTTCCTGAGCTTAAN 6FFFFFB?FBFFFBFFFFE@<@;FDFEC62,+;FEFFFFECFEFFFFFFFFFFFFFGFGFGGFGFGGGGFD>?7DGGGGGGGGGGGGGGFGGGGGGGGGGGEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGCCB8# NM:i:2 AS:i:602 RG:Z:1
not.rescuable 165 NC_001422.1 3017 0 * = 3017 0 GGTTTCGCAGCACCTGCACGGCGCCGTTTCTAGGACGTACATACCGGACAAGAGTATATCCGCGATGGACGGGTGATCCATGTAGGAAGTCTAGT GGGGGFFFFFFFFEEB@F29?3>EGGGGG8#GGGGGFAFEFFFFB=6>EEEGGGGBGGGGGGG;<;7EF;@EEC9CEFF??6?EF<@EDFFF### RG:Z:1
empty 77 * 0 0 * * 0 0 * RG:Z:1
empty 141 * 0 0 * * 0 0 * RG:Z:1
empty 77 * 0 0 * * 0 0 * * RG:Z:1
empty 141 * 0 0 * * 0 0 * * RG:Z:1
2 changes: 1 addition & 1 deletion tests/phix.se.sam
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -45,4 +45,4 @@ unmappable.41 4 * 0 0 * * 0 0 TGCATGCAGTACGTGCTCGCCGTAGGTTGTGGCTCATTCCTGCAACACGT
rescuable.42 16 NC_001422.1 3017 60 116=1X183=1X * 0 0 CTCCGGTGGCATTCAAGGTGATGTGCTTGCTACCGATAACAATACTGTAGGCATGGGTGATGCTGGTATTAAATCTGCCATTCAAGGCTCTAATGTTCCTAACCCTGATGAGGCCGTCCCTAGTTTTGTTTCTGGTGCTATGGCTAAAGCTGGTAAAGGACTTCTTGAAGGTACGTTGCAGGCTGGCACTTCTGCCGTTTCTGATAAGTTGCTTGATTTGGTTGGACTTGGTGGCAAGTCTGCCGCTGATAAAGGAAAGGATACTCGTGATTATCTTGCTGCTGCATTTCCTGAGCTTAAN 6FFFFFB?FBFFFBFFFFE@<@;FDFEC62,+;FEFFFFECFEFFFFFFFFFFFFFGFGFGGFGFGGGGFD>?7DGGGGGGGGGGGGGGFGGGGGGGGGGGEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGCCB8# NM:i:2 AS:i:602 RG:Z:1
rescuable.43 4 * 0 0 * * 0 0 ATAAGATCAGAAAATACAGCAGCAAAATAAACACGAGTATACTTTACTTTATCAGAGGCAAACTTACCACAAAGTACAACAAAATAAAGCAACTTATCAGAAACGACAGAAGTGCAAGCCAGCAAAGTACATTCAAGAAGACCTTAACCAACTTTAGCCAAAGCAACAGAAACAAAACTAAGGACAGCCTAATCAAGGTTAGGAAAATTAAAGCCATGAAAGGCAAATTTAATACAAGCAACACCAATGCATACAATATTATTATAGGTAACAAGAACATAACCTAGAATACCACAGGAG 8BCCGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGEGGGGGGGGGGGFGGGGGGGGGGGGGGD7?>DFGGGGFGFGGFGFGFFFFFFFFFFFFFEFCEFFFFEF;+,26CEFDF;@<@EFFFFBFFFBF?BFFFFF6 RG:Z:1
not.rescuable 16 NC_001422.1 3017 60 116=1X183=1X * 0 0 CTCCGGTGGCATTCAAGGTGATGTGCTTGCTACCGATAACAATACTGTAGGCATGGGTGATGCTGGTATTAAATCTGCCATTCAAGGCTCTAATGTTCCTAACCCTGATGAGGCCGTCCCTAGTTTTGTTTCTGGTGCTATGGCTAAAGCTGGTAAAGGACTTCTTGAAGGTACGTTGCAGGCTGGCACTTCTGCCGTTTCTGATAAGTTGCTTGATTTGGTTGGACTTGGTGGCAAGTCTGCCGCTGATAAAGGAAAGGATACTCGTGATTATCTTGCTGCTGCATTTCCTGAGCTTAAN 6FFFFFB?FBFFFBFFFFE@<@;FDFEC62,+;FEFFFFECFEFFFFFFFFFFFFFGFGFGGFGFGGGGFD>?7DGGGGGGGGGGGGGGFGGGGGGGGGGGEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGCCB8# NM:i:2 AS:i:602 RG:Z:1
empty 4 * 0 0 * * 0 0 * RG:Z:1
empty 4 * 0 0 * * 0 0 * * RG:Z:1
21 changes: 21 additions & 0 deletions tests/run.sh
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -69,3 +69,24 @@ if ${strobealign} --create-index tests/phix.fasta > /dev/null 2> /dev/null; then
# --details output is proper SAM
${strobealign} --details tests/phix.fasta tests/phix.1.fastq tests/phix.2.fastq 2> /dev/null | samtools view -o /dev/null
${strobealign} --details tests/phix.fasta tests/phix.1.fastq 2> /dev/null | samtools view -o /dev/null

# Secondary alignments

# No secondary alignments on phix
${strobealign} tests/phix.fasta tests/phix.1.fastq | grep -v '^@PG' > no-secondary.sam
${strobealign} -N 5 tests/phix.fasta tests/phix.1.fastq | grep -v '^@PG' > with-secondary.sam
test $(samtools view -f 0x100 -c with-secondary.sam) -eq 0
rm no-secondary.sam with-secondary.sam

# Secondary alignments for repeated phiX
cp tests/phix.fasta repeated-phix.fasta
echo ">repeated_NC_001422" >> repeated-phix.fasta
sed 1d tests/phix.fasta >> repeated-phix.fasta
${strobealign} repeated-phix.fasta tests/phix.1.fastq | grep -v '^@PG' > no-secondary.sam
${strobealign} -N 5 repeated-phix.fasta tests/phix.1.fastq | grep -v '^@PG' > with-secondary.sam
test $(samtools view -f 0x100 -c with-secondary.sam) -gt 0

# Removing secondary alignments gives same result as not producing them in the first place
samtools view -h --no-PG -F 0x100 with-secondary.sam > with-secondary-only-primary.sam
diff no-secondary.sam with-secondary-only-primary.sam
rm no-secondary.sam with-secondary.sam with-secondary-only-primary.sam repeated-phix.fasta