1
$\begingroup$

I need to identify all loci in the short-read sequence at which the number of microsatellite repeats (i.e. number of copies of "AA," "GTC," etc.) differ from the reference genome, as well as the locations at which the number of repeats is the same in both the reference genome and the short-read sequence.

I used the Burrows-Wheeler Aligner (BWA mem) to map a high-coverage short-read sequence (obtained from the NCBI's short-read sequence archive) to a reference genome. The output is in .sam format. I have also used a separate program to identify the loci in the reference genome at which microsatellites occur.

I would like to identify all the loci in the short-read sequence at which the microsatellite length and loci differ from the reference genome. Does anyone know any tools or packages I could use to read a .sam/.bam file of a short-read sequence mapped to a reference genome and identify specific loci at which the short-read sequence differs from the reference genome? I am using RStudio and have access to my university's supercomputer clusters.

For info on microsatellites, see here: https://en.wikipedia.org/wiki/Microsatellite#:~:text=A%20microsatellite%20is%20a%20tract,locations%20within%20an%20organism's%20genome.

My end goal is to determine the loci at which microsatellite alleles tend to differ within a species, as well as the loci at which microsatellite alleles tend to be identical between all members of the same species.

$\endgroup$
1
  • $\begingroup$ Can you please provide more information about your specific problem? This will greatly help in working through a solution. There is presumably some end-goal that this STR detection is used for. No algorithm will detect all variants, most will produce some amount of false negative and false positive results. $\endgroup$
    – gringer
    Jul 3 '20 at 2:33
3
$\begingroup$

I haven't done any whole-genome STR analysis from NGS data myself, but are aware of others that have used lobSTR for this. There's also a recent paper [here] that compares a few different STR analysis packages (i.e. RepeatSeq, LobSTR, HipSTR, GangSTR). Here's the concluding paragraph:

In conclusion, all these tools are built to genotype STRs but have different strengths and weaknesses. Based on our analysis there is no clear overall winner. RepeatSeq and HipSTR are the best when considering genotyping error rate even with low coverage. On the other hand, GangSTR has an advantage because it is the only tool among them that can call alleles longer than the read length but shows higher error rate unless looking at only the enclosed class of reads, which in turn would lose the GangSTR's advantage of picking up long genotypes. In addition, GangSTR is the newest tool and so comes with reference files for different reference builds that are periodically updated according to the tool's webpage. The correct choice of a tool and the subsequent filtering depends on the aim of the analysis, and might be influenced by available hardware resources and time limit for running tools.

In other words, without more information about the specific problem you're trying to solve, it's not clear what the best tool to use is (this is a common issue with Bioinformatics problems).

$\endgroup$
3
  • 1
    $\begingroup$ ah I see - thank you so much for linking the paper! I will give it a read. My specific problem is that I need to identify all loci in the short-read sequence at which the number of microsatellite repeats (i.e. number of copies of "AA," "GTC," etc.) differ from the reference genome, as well as the locations at which the number of repeats is the same in both the reference genome and the short-read sequence $\endgroup$ Jul 2 '20 at 22:21
  • $\begingroup$ Honestly, if anyone knows of a software that works on R and detects indels using BWA-generated files, that would work for my purposes $\endgroup$ Jul 2 '20 at 22:45
  • $\begingroup$ Yep its a cool answer, I like it too. Your comment above is a separate question. My personal thoughts are you'd hit a RAM bottleneck pretty quickly via R and you should think about a more complex data pipeline. $\endgroup$
    – M__
    Jul 3 '20 at 1:50

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.