<p>I'm trying to study fetal fraction by using a seqff script from <a href="https://obgyn.onlinelibrary.wiley.com/doi/abs/10.1002/pd.4615" rel="nofollow noreferrer">here</a>. The instruction says that I need a headless sam file to work so after alignment with bwa, I process my sam file by using samtools:</p>
<pre><code>samtools view -S -q 20 -F 4 target.sam &gt; treated.sam
</code></pre>
<p>This process also remove the @SG header since there is no reference file. I also use grep to remove header of the target.sam got after alignment for comparison. However, when I compare, I noticed that ENet calculation on the treated.sam returned with NA results, so half of the method to calculate fetal fraction of this script was gone. I have tested it in multiple cases and it seems that without using samtools, the script worked well but with samtools Enet calculaton not worked. Can someone tell me what is wrong with the script? I'm using R 4.0 and all modification I made to the script is on the args since it cannot recognize args I input with anymore.</p>
<p>p/s: this is few lines of headless target.sam made by using grep</p>
<pre><code>W2003120038L1C001R0050000009    0   chr10   6490952 0   50M *   0   0   GAGCAGAAAAACTGGAAACTCTAAAAAGCAGAGTGCCTCTCCTCCTCCAA  IBGIBAB&lt;ABBCA9IB@BHDIA=ABBB8IB9A;DGII@I8IIAII&lt;IDCA  XT:A:R  NM:i:0  X0:i:151    XM:i:0  XO:i:0  XG:i:0  MD:Z:50
W2003120038L1C001R0050000037    0   chr5    172722817   37  50M *   0   0   GTCCTAAGACAAAGACCTCAGTGCTTATAATGTAAGTGGTCACTGGGCTG  ICIICBB:BIBBBIAICDIAIAIICDBDAABIDBBICEIAIAFDH?IIAI  XT:A:U  NM:i:1  X0:i:1  X1:i:0  XM:i:1  XO:i:0  XG:i:0  MD:Z:7T42
W2003120038L1C001R0050000062    0   chr5    84724587    37  50M *   0   0   TAAGAATTGTGCATAACCAGTAGCCCATGTTATCTTGCAGAAGCATTTGT  ABBI@BDDIDIIADABHEAIC@IIIBBAID@ADIDCA:BC=BIIAAAAIC  XT:A:U  NM:i:0  X0:i:1  X1:i:0  XM:i:0  XO:i:0  XG:i:0  MD:Z:50
W2003120038L1C001R0050000067    16  chr5    45648067    37  50M *   0   0   TTTATTTAAATACAAGCCAATTAAATGATTTTATTTTACATTTTCTAACC  AABAAABCAAABH@CFHIDDB@BCDAID@BBACBB@BBIDB@2BIBC@II  XT:A:U  NM:i:0  X0:i:1  X1:i:0  XM:i:0  XO:i:0  XG:i:0  MD:Z:50
W2003120038L1C001R0050000184    16  chr4    103629659   37  50M *   0   0   CAGTGCTAAGAGGGAAGTTTATAGCAACAGATACTTACATCCAAAAAGTA  FDFBIIADAFDICICCIAABDBAFIADIDFDABIBA@I@CI3?DCCDIBC  XT:A:U  NM:i:0  X0:i:1  X1:i:0  XM:i:0  XO:i:0  XG:i:0  MD:Z:50
W2003120038L1C001R0050000211    0   chr3    157734199   37  50M *   0   0   TCTGAAATTGATCTTTATCCCACATCCACGTATCCTAAAATAATAAAACA  AFDBBBB*CIBDFACABDIIBB&lt;1CI=BFIDBAFEABBAADBBABAABIB  XT:A:U  NM:i:2  X0:i:1  X1:i:0  XM:i:2  XO:i:0  XG:i:0  MD:Z:7G15T26
</code></pre>
<p>And this is from treated.sam, the sam file be processed by samtools</p>
<pre><code>W2003120038L1C001R0050000037    0   chr5    172722817   37  50M *   0   0   GTCCTAAGACAAAGACCTCAGTGCTTATAATGTAAGTGGTCACTGGGCTG  ICIICBB:BIBBBIAICDIAIAIICDBDAABIDBBICEIAIAFDH?IIAI  XT:A:U  NM:i:1  X0:i:1  X1:i:0  XM:i:1  XO:i:0  XG:i:0  MD:Z:7T42
W2003120038L1C001R0050000062    0   chr5    84724587    37  50M *   0   0   TAAGAATTGTGCATAACCAGTAGCCCATGTTATCTTGCAGAAGCATTTGT  ABBI@BDDIDIIADABHEAIC@IIIBBAID@ADIDCA:BC=BIIAAAAIC  XT:A:U  NM:i:0  X0:i:1  X1:i:0  XM:i:0  XO:i:0  XG:i:0  MD:Z:50
W2003120038L1C001R0050000067    16  chr5    45648067    37  50M *   0   0   TTTATTTAAATACAAGCCAATTAAATGATTTTATTTTACATTTTCTAACC  AABAAABCAAABH@CFHIDDB@BCDAID@BBACBB@BBIDB@2BIBC@II  XT:A:U  NM:i:0  X0:i:1  X1:i:0  XM:i:0  XO:i:0  XG:i:0  MD:Z:50
W2003120038L1C001R0050000184    16  chr4    103629659   37  50M *   0   0   CAGTGCTAAGAGGGAAGTTTATAGCAACAGATACTTACATCCAAAAAGTA  FDFBIIADAFDICICCIAABDBAFIADIDFDABIBA@I@CI3?DCCDIBC  XT:A:U  NM:i:0  X0:i:1  X1:i:0  XM:i:0  XO:i:0  XG:i:0  MD:Z:50
W2003120038L1C001R0050000211    0   chr3    157734199   37  50M *   0   0   TCTGAAATTGATCTTTATCCCACATCCACGTATCCTAAAATAATAAAACA  AFDBBBB*CIBDFACABDIIBB&lt;1CI=BFIDBAFEABBAADBBABAABIB  XT:A:U  NM:i:2  X0:i:1  X1:i:0  XM:i:2  XO:i:0  XG:i:0  MD:Z:7G15T26
W2003120038L1C001R0050000212    0   chr4    157022243   37  50M *   0   0   AAAAATAAAAAAGCATTCTGCCCAGATCTCTCAACATTATGTTATTTCTT  4&lt;AA?D2,AAA@IF@DDIDIIDIA2AD=DBAF@BIADA&lt;DICD@DADFD@  XT:A:U  NM:i:0  X0:i:1  X1:i:0  XM:i:0  XO:i:0  XG:i:0  MD:Z:50
W2003120038L1C001R0050000258    0   chr13   106093727   37  50M *   0   0   AATGGTTACTTTATGTTTATATTTCAACAGGGAATATTTGAAATATGGCA  97?IIDA*FDAD&lt;@H@AD&lt;@&lt;?AAG8??,I;I@:DADDAG7@@7(;IH6?  XT:A:U  NM:i:1  X0:i:1  X1:i:0  XM:i:1  XO:i:0  XG:i:0  MD:Z:7C42
</code></pre>

CC BY-SA 4.0
2021-04-22T07:35:12.350Z
@gringer I have added few lines of both file on my post. Also I have to perform quality filter to filter out bad alignments results
khanhlpbao
CC BY-SA 4.0
2021-04-22T11:26:24.447Z
If the grep one works, and the non-grep one doesn't, my first step would be to exclude any confounding factors (like quality filtering) to try to get them both to look the same (which they don't, at the moment). Is there any improvement if you don't include the quality filter (either in the files looking similar, or the script working better)?
gringer
CC BY-SA 4.0
2021-04-22T23:11:14.637Z
If you think there is a formatting problem, you might need to put small subsets of the samtools file through until you can identify the problematic line.
swbarnes2
CC BY-SA 4.0
2021-04-23T04:01:33.283Z
@swbarnes2 thanks, do you know any command or program that can check for the different? Because the each sam file is more than 4GB, finding what is different by compare with text editor is impossible to say at least
khanhlpbao
CC BY-SA 4.0
2021-04-23T17:02:40.410Z
One you figure out the problem in one file, you should be able to fix it for the others.
swbarnes2
<r><samtools>
seqff script got trouble after samtools treatment
# seqff script got trouble after samtools treatment

(Asked by: username_0 on 2021-04-22)

I'm trying to study fetal fraction by using a seqff script from here. The instruction says that I need a headless sam file to work so after alignment with bwa, I process my sam file by using samtools:
samtools view -S -q 20 -F 4 target.sam > treated.sam

This process also remove the @SG header since there is no reference file. I also use grep to remove header of the target.sam got after alignment for comparison. However, when I compare, I noticed that ENet calculation on the treated.sam returned with NA results, so half of the method to calculate fetal fraction of this script was gone. I have tested it in multiple cases and it seems that without using samtools, the script worked well but with samtools Enet calculaton not worked. Can someone tell me what is wrong with the script? I'm using R 4.0 and all modification I made to the script is on the args since it cannot recognize args I input with anymore.
p/s: this is few lines of headless target.sam made by using grep
W2003120038L1C001R0050000009    0   chr10   6490952 0   50M *   0   0   GAGCAGAAAAACTGGAAACTCTAAAAAGCAGAGTGCCTCTCCTCCTCCAA  IBGIBAB
And this is from treated.sam, the sam file be processed by samtools
W2003120038L1C001R0050000037    0   chr5    172722817   37  50M *   0   0   GTCCTAAGACAAAGACCTCAGTGCTTATAATGTAAGTGGTCACTGGGCTG  ICIICBB:BIBBBIAICDIAIAIICDBDAABIDBBICEIAIAFDH?IIAI  XT:A:U  NM:i:1  X0:i:1  X1:i:0  XM:i:1  XO:i:0  XG:i:0  MD:Z:7T42
W2003120038L1C001R0050000062    0   chr5    84724587    37  50M *   0   0   TAAGAATTGTGCATAACCAGTAGCCCATGTTATCTTGCAGAAGCATTTGT  ABBI@BDDIDIIADABHEAIC@IIIBBAID@ADIDCA:BC=BIIAAAAIC  XT:A:U  NM:i:0  X0:i:1  X1:i:0  XM:i:0  XO:i:0  XG:i:0  MD:Z:50
W2003120038L1C001R0050000067    16  chr5    45648067    37  50M *   0   0   TTTATTTAAATACAAGCCAATTAAATGATTTTATTTTACATTTTCTAACC  AABAAABCAAABH@CFHIDDB@BCDAID@BBACBB@BBIDB@2BIBC@II  XT:A:U  NM:i:0  X0:i:1  X1:i:0  XM:i:0  XO:i:0  XG:i:0  MD:Z:50
W2003120038L1C001R0050000184    16  chr4    103629659   37  50M *   0   0   CAGTGCTAAGAGGGAAGTTTATAGCAACAGATACTTACATCCAAAAAGTA  FDFBIIADAFDICICCIAABDBAFIADIDFDABIBA@I@CI3?DCCDIBC  XT:A:U  NM:i:0  X0:i:1  X1:i:0  XM:i:0  XO:i:0  XG:i:0  MD:Z:50
W2003120038L1C001R0050000211    0   chr3    157734199   37  50M *   0   0   TCTGAAATTGATCTTTATCCCACATCCACGTATCCTAAAATAATAAAACA  AFDBBBB*CIBDFACABDIIBB<1CI=BFIDBAFEABBAADBBABAABIB  XT:A:U  NM:i:2  X0:i:1  X1:i:0  XM:i:2  XO:i:0  XG:i:0  MD:Z:7G15T26
W2003120038L1C001R0050000212    0   chr4    157022243   37  50M *   0   0   AAAAATAAAAAAGCATTCTGCCCAGATCTCTCAACATTATGTTATTTCTT  4

**Tags:** r, samtools
> **Comments:**
> **username_0**: @gringer I have added few lines of both file on my post. Also I have to perform quality filter to filter out bad alignments results  
> **username_1**: If the grep one works, and the non-grep one doesn't, my first step would be to exclude any confounding factors (like quality filtering) to try to get them both to look the same (which they don't, at the moment). Is there any improvement if you don't include the quality filter (either in the files looking similar, or the script working better)?  
> **username_2**: If you think there is a formatting problem, you might need to put small subsets of the samtools file through until you can identify the problematic line.  
> **username_0**: @swbarnes2 thanks, do you know any command or program that can check for the different? Because the each sam file is more than 4GB, finding what is different by compare with text editor is impossible to say at least  
> **username_2**: One you figure out the problem in one file, you should be able to fix it for the others.  

---

## 0 answers

