Questions tagged [samtools]

Questions specific to interacting with and post-processing sequence alignments using the SAMtools package. For questions specifically about formats SAM, BAM or CRAM use tags sam or bam,

Filter by
Sorted by
Tagged with
1
vote
0answers
37 views

Conversion of SAM to BAM files

I am very new to micro RNA analysis. I have been using H. sapiens, GRCh38 + major index as given in the Bowtie Website to align with my trimmed FASTQ file . The command I am using to make very ...
2
votes
1answer
21 views

Python package or CLI tool outputting mutations of sequence with respect to reference

I'm new to bioinformatics, and I'm looking for a tool that takes two FASTA files, one containing a reference, let's call it A.fasta, the other another sequence, let'...
0
votes
0answers
29 views

seqff script got trouble after samtools treatment

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: ...
1
vote
1answer
40 views

Inflate operation failed: progress temporarily not possible, or in() / out() returned an error

I download a BAM file from ftp-trace.ncbi.nlm.nih.gov:/giab/ftp/data/NA12878/10Xgenomics_ChromiumGenome_LongRanger2.0_06202016/NA12878_GRCh38.bam (171G), and I am trying to view the comparison result ...
0
votes
0answers
46 views

bcftools mpileup double counting for overlapping reads

I have a dataset of paired-end reads, which can overlap. So when i run bcftools mpileup this overlap gives a double count in the depths so we see DP = 2 even though it should be DP = 1. The question ...
0
votes
1answer
28 views

How can I run the latest version of samtools on Slurm?

I am trying to run the following script on slurm to extract information from .sam files. While everything is perfect on my local machine with the samtools version 1.11, it does not work on slurm. I ...
0
votes
1answer
53 views

How to perform bowtie2 analysis with slurm?

I am trying to run my alignment script that works locally, using SBATCH. Official manual for bowtie2 says I can use ...
1
vote
1answer
93 views

How can I extract information from .sam files?

I have 10 .sam files after my bowtie2 alignment on ten single-pair sequences. I would like to build a graph based on that output ...
1
vote
1answer
245 views

What is “unmapped read segments” in the output of samtools idxstats?

samtools idxstats produces a four column output (see here) ...
1
vote
1answer
25 views

How to extract reads with INDELs > a given size?

I'm trying to modify this https://www.biostars.org/p/253774/ To get reads with deletions > 20bp I think this gives reads with exactly 20bp dels: ...
2
votes
1answer
152 views

How to generate a consensus sequence from a multi-reference BAM file?

I am trying to generate a consensus sequence from a BAM file that was generated by mapping reads to a reference FASTA containing multiple sequences. Usually, I generate consensus sequences from BAM ...
0
votes
1answer
74 views

Reading fasta sequence regions with faidx

Im coding in c++ and reading in different reference genomes to examine regions across the chromosomes a few hundred basepairs at a time. To do it in c++ i use the library htslib and the command ...
0
votes
1answer
26 views

How to identify to each scaffold a read belongs to, inside a .sam file?

I have a fasta file assembly and combining it with the raw reads we produced a .bam file which I converted to .sam . The .sam information lines look like this: ...
1
vote
2answers
173 views

Trying to create a .bam file without the need for a .sam file

I'm trying to use the code specified in this link to create a .bam file without the need for a .sam file. Here is the code I'm using: ...
1
vote
2answers
87 views

Mapping statistics from bam file using bbtools and sambamba

Below are the statistics for RNA-seq mapped and unmapped paired-end reads to rice genome using reformat.sh from bbtools on bam files. It gives 77% mapped and 5% unmapped, what about the remaining 18% ...
0
votes
0answers
49 views

summarising read group information from a .bam file

I have merged together 2 different .bam files in order to simulate sample contamination. So the reads can come from one of two samples, as shown by the read group info: ...
0
votes
1answer
170 views

NO_COOR reads not in a single block at the end 0 -1

I am doing an RNA seq analysis with mouse samples. I had my reads in several lanes, so I concatenated those fastq reads. Then I've mapped my reads with Hisat2 against the new version of mouse "...
0
votes
0answers
9 views

How do I tell IsoMut2py what to input into mpileup?

Just a warning, this question is extremely specific to the indel-detection software IsoMut2py. The paper for isomut can be found here: https://bmcbioinformatics.biomedcentral.com/articles/10.1186/...
0
votes
2answers
47 views

Basic question about finding mRNA sequence in transcriptome

I am a computer scientist just starting out in bioinformatics topics, and would appreciate any guidance that can be given here: I have an mRNA sequence- an isoform - whose length is about 4000 base ...
0
votes
1answer
1k views

Samtools Index: Chromosome Blocks not Continuous

I am working with short-read whole-genome sequences from the NCBI's SRA. I have aligned and sorted all of my short-read sequences and am attempting to index each sequence into .bai format using ...
0
votes
0answers
31 views

How do I find all supplementary mappings of a chimeric long (ONT) read mapped with NGMLR?

Unlike this guy How do I find split reads? I have done a lot of reading before asking this question. I'm aware that the FLAG will tell me if a read is supplementary and the CIGAR has information ...
2
votes
1answer
207 views

How to extract unmatched reads using bwa and samtools?

I have a single read (NOT paired) that I need to pass through the workflow described in Beauclair et al. paper (free version here https://rnajournal.cshlp.org/content/24/10/1285.long) for identifying ...
1
vote
1answer
58 views

How can I get unmatched reads for defective genomes analysis using bwa and samtools?

I am trying to follow the Beauclair et al. paper (free version here https://rnajournal.cshlp.org/content/24/10/1285.long) for identifying defective genomes using their DI-tector program. According to ...
2
votes
4answers
100 views

How to make work programs from the $PATH?

I am trying to analyze my RNAseq reads for defective genomes and I use this program (http://www.di-tector.cyame.eu/) that is suggested by Beauclair et al (https://pubmed.ncbi.nlm.nih.gov/30012569/) ...
3
votes
1answer
128 views

Is there a tool that can perform a read-group-aware mpileup from a single file?

I would like to perform a samtools mpileup from a single file that contains thousands of read groups with different SM tags. I could split the bam by read group ...
2
votes
1answer
30 views

Phasing problem using “samtools phase”

so I'm having a problem using samtools phase. I'm trying to phase my .bam archives that where generated by aligning my samples with my reference genome. Following the documentation I tried the command:...
0
votes
2answers
247 views

Accessing .bam/.cram files from AWS S3?

What's the best/easiest way for accessing .bam/.cram files from S3? I have .bam/.bai .cram/.crai and .bed + .gff files that I want to make available easily to others so they can browse in IGV and ...
1
vote
2answers
1k views

merging two/or more bed file into one bed file

I am trying to merge two bed files (more in future) to one. my bed files are something like : . I need to merge them in a way to have the shared chromosome location. Is there a way to do that ?
1
vote
1answer
162 views

Cuffmerge: EOF marker is absent. The input is probably truncated

I need to merge my all transcripts.gtf from cufflinks output follow this command line : cuffmerge -o merged_gtf_output -p 15 -s ref.fasta -g anot.gtf assembly.txt ...
1
vote
2answers
388 views

Marking optical or PCR duplicates with picard vs. samtools flagstat

I am trying to identify the best library prep method for noninvasive prenatal test samples, to be sequenced with Illumina Novaseq. One metric that I am evaluating is the number/ percentage of PCR or ...
1
vote
0answers
19 views

Indel quality from mpileup

Could someone please confirm that I understand this correctly? I have the res and qual columns from mpileup and I would like to match them to get the qual per base. It seems that the indel initiation (...
2
votes
1answer
70 views

Why does samtools mpileup sometimes include ref bases (other than ',' or '.')?

This is my first post here. I can think of no way of giving an easily reproducible example, as per the stack ethos. SO apologies in advance and any feedback on question format appreciated. I have ...
1
vote
1answer
320 views

How to filter a SAM file by a bed file?

I have a large BAM file and I want to filter it by just what's defined in a bed file, and write it to a new file. Is it possible to do this with samtools, and if so, how? Is there a better tool I ...
1
vote
1answer
61 views

samtools / bamUtil | Meaning of <wildcard> as Reference Name

I have been working with several tools on bam files recently and I am not sure how I should interpret some of the outputs. samtools idxstats mapped.sorted.bam <...
1
vote
1answer
144 views

How to obtain file offsets corresponding to range query in block-gzipped file

How can we retrieve the file offsets of the zipped-chunks in a block-gzipped file that contain records for a region/window of interest (e.g. Chr01:100-Chr02:145)? It seems like the CLI of ...
3
votes
1answer
234 views

How to generate two fasta files from samtools phase output?

samtools phase is an easy-to-use tool to phase variants from a bam file and a reference. For a sequence that has two haplotypes, like this: ...
0
votes
1answer
124 views

Single-cell sequencing dataset has too many barcodes

I am analyzing a single-cell sequencing dataset from the website 10xgenomics, with 2000 cells. It is a BAM file and I am trying to obtain the individual cells per sample. I used the command ...
0
votes
1answer
47 views

Obtaining HGDP project data in fasta format

I need to obtain sample data from modern humans in fasta format. I just need some megabytes of data from every individual. I actually use a script that obtains the cram file from here (ftp.1000genomes....
3
votes
1answer
197 views

pysam or piping samtools view to a python script

Is it faster to use the PySam package to run a python script on a bam for read in samfile.fetch('chr1', 100, 120): print read compared to using a pipe and ...
0
votes
1answer
124 views

How to convert a Pileup file to VCF format with Hg19 alignment

I would like to know how to convert a pileup file for an adna y chromosome to vcf with Hg19 alignment.
3
votes
2answers
253 views

What is a samtools mpileup reference skip?

The samtools documentation for mpileup states: At this column, a dot stands for a match to the reference base on the forward strand, a comma for a match on the reverse strand, a '>' or '<' for a ...
1
vote
2answers
366 views

bedtools intersect on very large .bam file - -sorted confusion

I need to do a bedtools intersect operation on a very large .bam file. When I use the standard bedtools intersect operation, the process consumes all the memory on my system. From the bedtools ...
1
vote
4answers
70 views

Find indels between two short sequences

I have two sequences, say AAAGCTCGAGG and AAAGCGAGG. I need a convenient tool which shows me insertions or deletions between these, i.e. in this case something like ...
2
votes
1answer
161 views

How to take into account alternative bwa mem mapping when computing coverage

when mapping short reads with bwa mem, if a read has alternative mapping positions they are reported by bwa mem in the X0 and XA tags. Now, let's say I want to compute the coverage of my bam file. ...
5
votes
2answers
234 views

samtools mpileup skipping read

I run the command: samtools mpileup -O -s -q20-B -Q20 -f hg19.fa -r chr1:569929-569931 myFile.bam and get: ...
3
votes
1answer
345 views

How to get a MSA fasta from BAM/SAM?

I am working with paired end NGS miseq data from a viral genome with multiple timepoints. I have filtered and trimmed this data for quality and adapter sequences. I have then merged the filtered ...
6
votes
1answer
1k views

Convert BAM to properly paired FASTQ files

I thought I had figured this one out. But apparently not. I need to convert a BAM file of paired-end alignments to two FASTQ files of paired reads to realign them, with a twist: I only want reads ...
3
votes
2answers
2k views

Error while calling bcftools mpileup - Failed to open -: unknown file type

I have sequenced a bacterial genome with a GridIon from ONT. Basically what I want to check is whether or not trimming 50 bps at the beginning of the reads will improve alignment against the reference ...
2
votes
3answers
752 views

How to reverse complement the DNA sequences for given inverse/reverse coordinates?

I have the series of coordinates in id.txt file, whose coordinates sequences are in genome.fasta file. The coordinates of id.txt ...
3
votes
2answers
242 views

Frequency of specific viral sequence in .BAM or .fastq

I was wondering if anyone had any experience 'counting' the frequency of a particular viral sequence in an individual's fastq sequence. So basically, counting the number of occurrences of a particular ...