Hello,
I apologize if this forum is not the preferred place for this question, but are there standardized methods to infer the amplicons given a Bam/Sam file? That is, given a bam/sam file as input, how can I infer the starting and ending positions of the amplicons? I am dealing with viral sequence data which can be highly variable, so start and end positions may vary and partitioning may be difficult.
Edit: For further clarification, suppose amplicon-based sequencing was done over 3 separate PCR experiments. The experiments were designed in such a way so that the first overlaps with the second, and the second with the third. The output of each run is aggregated into a single FASTA file and then mapped to some reference which produces a bam/sam file. Now, it is my job to infer the sepratate original experiments. That is, I must take the bam file and partition it to produce 3 sub-datasets, each corresponding to the original PCR windows. Given that start/end positions within an amplicon aren't the same, what is the best/standard approach?
0 answers
No answers yet.
Log in to answer this question.
-
Issue softclipping reads when they belong and don't belong to a common amplicon
written by ManuelDB
-
retrieve amplicon sequences from fna file
written by Lila M
-
BAM header edit with a list of amplicons
written by ltalignani •
-
remove primer sequences from BAM file
written by J.F.Jiang
-
extract region bam
written by J.F.Jiang
-
VarDict Amplicon-Aware Calling with barcoded seq data
written by jsneaththompson
-
Need to parse pooled bam file based on sequence
written by djvanderweele •
-
Color / mark specific sequence bases differently in IGV
written by Matt Shirley
-
Query Ucsc With Known Sequence To Get Chr Position
written by mathew.nightingale •
-
How Are Ensembl'S Gene Start And End Coordinates Determined?
written by Dr.Bunsen •
amplicon ? do you mean "segment"= a pair of short reads ?
The reads are amplicon-based, rather than shotgun-based. I have a tool I developed that does haplotype reconstruction of viral sequences given amplicon-reads. However the tool depends on knowing the start and end positions of each amplicon experiment (i.e. where they should map to reference). But I would like the final product (for usability reasons) to just take as input a bam/sam file, rather than a FASTA file and starting/ending positions for each amplicon.
So you want to know where the reads in the sam/bam file map to the reference amplicon?
Well, slightly more than that. I would like to partition the data into what would have been each individual amplicon experiment.
Sorry. I am still confused on what you are trying to do here. What did you map your reads to exactly?
I have edited the original question to try and clarify.
So you mapped to the 3 amplicon sequences?
Yes. Is my assumption that start/end positions within an amplicon after alignment will be "jagged" in a sense correct? I guess you could trivially identify the regions after mapping if "-" characters are used for low coverage-regions in the beginning and end, but what about proper covered regions?
If you mapped to the 3 amplicon sequences, why not just partition your sam/bam file by which amplicon the read mapped to?
Is it reasonable to expect that each experiment will be aggregated into a single BAM file? Or is it ok to have an individual BAM file for each experiment. That would avoid this entire problem, but at more work for the user.
I'm also confused by your question. But I have a basic question: do the PCR windows overlap? You mention doing 3 separate experiments with overlap, and combining the results together (but not in such a way that the results were tagged somehow so they could be later de-convoluted after mixing). If the PCR windows (i.e. the amplicons) overlap, I'm not sure how you would be able to easily deconvolute the results. Depending on sequence depth, you might be able to see different heights for different amplicons, and infer something that way. Do I understand correctly that an experiment consists of PCR amplification of viral sequence - followed by by high-throughput sequencing of the amplicon (1 or many?), and this was done three times for amplicons with some degree of overlap?
Sequence depth/coverage is the way I had currently planned on doing it. Overlap positions should have very large peaks in depth/coverage. By finding these quasi-global maxima we can identify a starting/ending area for overlapping amplicons. The only thing left is starting of first amplicon and ending of last. That can be done by finding the local maxima as well (simple march forward/backward until good enough coverage).