55

Let’s start with what they have in common: All three formats store sequence data, and sequence metadata. Furthermore, all three formats are text-based. However, beyond that all three formats are different and serve different purposes. Let’s start with the simplest format: FASTA FASTA stores a variable number of sequence records, and for each record it ...


41

Genomes are commonly stored as either fasta files (.fa) or twoBit (.2bit) files. Fasta files store the entire sequence as text and are thus not particularly compressed. twoBit files store each nucleotide in two bits and contain additional metadata that indicates where there's regions containing N (unknown) bases. For more information, see the ...


25

The standard formats for storing sequence data are fasta and fastq. Fasta is used if you only need the raw sequence data, fastq is used if you want to store the sequence data along with the quality information from base calling. Each of these can be compressed using gzip or another standard compression algorithm. Typically we want to keep the quality ...


23

You're the second person I have ever seen using NCBI "chromosome names" (they're more like supercontig IDs). Normally I would point you to a resource providing mappings between chromosome names, but since no one has added NCBI names (yet, maybe I'll add them now) you're currently out of luck there. Anyway, the quickest way to do what you want is to samtools ...


20

The standard and the most common sequence format is FASTA for sure. You can compress it with a compressor. For the ~3GB human genome, gzip reduces the size to ~900MB, depending on the option in use. Another often used format is UCSC's 2-bit format. This format keeps each A/C/G/T with 2 bits. As I remember, it keeps non-A/C/G/T bases and lowercases in two ...


18

I'm not sure I'm doing it the best way, but here is an example where I read a compressed gzip fastq file and write the records in block gzip fastq: from Bio import SeqIO, bgzf # Used to convert the fastq stream into a file handle from io import StringIO from gzip import open as gzopen records = SeqIO.parse( # There is actually simpler (thanks @peterjc) ...


17

This is easy to check, you can download both specs in .tex format and do diff. Changes to the v4.2 compared to v4.1: Information field format: adding source and version as recommended fields. INFO field can have one value for each possible allele (code R). For all of the ##INFO, ##FORMAT, ##FILTER, and ##ALT metainformation, extra fields can be included ...


17

In a nutshell, FASTA file format is a DNA sequence format for specifying or representing DNA sequences and was first described by Pearson (Pearson,W.R. and Lipman,D.J. (1988) Improved tools for biological sequence comparison. Proc. Natl Acad. Sci. USA, 85, 2444–2448) FASTQ is another DNA sequence file format that extends the FASTA format with the ability ...


17

Whenever you want to save space (this can be a substantial savings). Until quite recently (samtools/htslib 1.7), only CRAM supported long CIGAR strings. If you need to guarantee that any random obscure downstream program will be able to handle it. Uptake of CRAM has been pretty slow. Java programs using htsjdk (e.g., picard, IGV and GATK) have only ...


14

You can do this easily with bioawk, which is a version of awk with added features facilitating bioinformatics: bioawk -c fastx '{print $name"\t0\t"length($seq)}' test.fa -c fastx tells the program that the data should be parsed as fasta or fastq format. This makes the $name and $seq variables available in the awk commands.


13

There are several things to consider when asking for "the most efficient" way to store data, it all depends on your use case. Do you just need ACGT, or are there also IUPAC codings for combinations? Do you need additional data (like quality values)? What kind of application are you using the data for (does it need to load all at once or in in chunks? Once or ...


13

Arithmetic on a zero-based coordinate system is less complicated than that on a one-based system, so it appears zero-based is often (not exclusively) used for binary data formats, like BAM or bigBed, or text formats like BED, where computers are used to more efficiently calculate lengths of or set operations on intervals. A more complete answer on the ...


11

This can be done in R very easily from an indexed .bam file. Given single-end file for sample1. library(GenomicAlignments) library(rtracklayer) ## read in BAM file (use readGAlignmentPairs for paired-end files) gr <- readGAlignments('sample1.bam') ## convert to coverages gr.cov <- coverage(gr) ## export as bigWig export.bw(gr.cov,'sample1.bigwig') ...


11

It's good practice to have your FASTA indexed, so you can leverage the .fai you are likely to already have. If not, you can just generate the index with samtools and use some awk to make your BED: samtools faidx $fasta awk 'BEGIN {FS="\t"}; {print $1 FS "0" FS $2}' $fasta.fai > $fasta.bed This will maintain tab separation but you can drop the BEGIN ...


11

Having done 23andme myself I can tell you that your variant file, which contains SNP genotypes, cannot be converted to a bam file. It does not contain the same information as a bam file. It may be helpful to familiarize yourself with those filetypes and the technologies used to obtain them. A SAM/BAM file contains alignments of reads obtained by sequencing. ...


9

The sample tag (i.e. SM) was a mandatory tag in the initial SAM spec (see the .pages file; you need a mac to open it). When transitioned to Latex, this requirement was mysteriously dropped. Picard is conforming to the initial spec. Anyway, the sample tag is important to quite a few tools. I would encourage you to add it.


8

I think the question is a bit ambiguous so please excuse this answer that's a bit redundant from the rest of the ones provided. As others have mentioned, if you want to store a full genome, FASTA and 2bit formats are appropriate. For some context, hg19 is about 900Mb compressed for the FASTA file and about 780Mb compressed for the 2bit file . hg19 is a ...


8

The first place to start is the GFF3 specification. This is the official word on what is and is not allowed in a GFF3 file. For example, users can define arbitrary attribute keys, so long as they do not begin with an uppercase letter (these are reserved for "official" use). But your question doesn't seem to be about what is allowed, but what is commonly ...


7

It is not yet standardized, but graph format has the potential for being the most space-efficient method for storing genomes. The idea is this: rather than store a genome as a linear string of sequenced nucleotides, genomes are stored as overlapping graphs, where sequence variants branch off from the reference genome, and then rejoin when the alignment ...


7

FASTA (officially) just stores the name of a sequence and the sequence, unofficially people also add comment fields after the name of the sequence. FASTQ was invented to store both sequence and associated quality values (e.g. from sequencing instruments). SAM was invented to store alignments of (small) sequences (e.g. generated from sequencing) with ...


7

You mention Biopython, which contains tests: https://github.com/biopython/biopython/tree/master/Tests. Some of the tests consist in reading files present in the folders listed in the above link. These files could be a starting point for a database of test files. Whenever one comes across a test case not covered with these files, one could construct a new ...


7

This is an early Solexa/Illumina sequencer format. The columns are the identifier to location on the flowcell. I believe the "." was the original placeholder for an unread base, which has been replaced by and "N" in current Illumina sequencing output. From on http://www.crg.eu/en/content/processing-and-analysis-illumina-sequencing-data. ...


7

This should work for you: for i in *.tar.gz do dir2="/destination/directory" tar xvzf $i -C $TMPDIR files=($TMPDIR/*) module load HISAT2/2.0.4-goolf-1.7.20; module load SAMtools/1.3.1-goolf-1.7.20; hisat2 -p 8 --dta --rna-strandness RF -x /genome_index -1 $files[0] -2 $files[1]| samtools view -Sb - > ${dir2}/$i.bam done (If this is the ...


6

You could adapt this awk one-liner. Note that it assumes that sequence IDs are not longer than 100 characters and that there is no description following the sequence ID on the header line. cat myseqs.fasta | awk '$0 ~ ">" {print c; c=0;printf substr($0,2,100) "\t0\t"; } $0 !~ ">" {c+=length($0);} END { print c; }' Otherwise, any Bio* library (Perl, ...


6

By default, a CRAM you create with samtools is lossless. It typically halves the input BAM in terms of file size. If you want to compress more, you can let samtools convert most read names to integers. You won't be able to tell optical duplicates from read names, but this is a minor concern. You can also drop useless tags depending on your mapper and the ...


6

Not as far as I am aware. The Ray assembler used to (and possibly still does) store the kmers as FASTA files where the header was the count of the sequence, which I thought was a pretty neat bastardisation of the FASTA file format. It looks like this format is also used by Jellyfish when reporting kmer frequencies by the dump command (but its default output ...


6

It looks like a weird way to represent data. You're right, it does not at all look like a BED. Digging in GEO made me find this information: Supplementary_files_format_and_content: Tab-delimited file reports one contact per row. Interacting HindIII fragments are represented in bed format (chr/start/stop) with bait (or upstream bait, in the case of ...


6

So binaryCIF is newer and is heavily inspired by MMTF. NGL, the most "loved" JS library (IMO) to show protein fetches proteins as MMTF, while Mol* uses bCIF —Mol* has the official repo of the bCIF format. Alex Rose, developed the former at the RCSB PDB, and the latter in collaboration with the PDBe. He also worked with Antony Bradley, the main ...


5

Not that I am aware. It is best to go with format specifications when coding. Also it may be good to look at the example files that come together with various tools performing file conversions and handling. E.g. Trimmomatic comes with few .fa files (http://www.usadellab.org/cms/?page=trimmomatic) Samtools with a toy.sam (https://github.com/lh3/samtools-...


5

I would consider the description there a bug. The filter is actually the strand, strand is the frame, group is the attribute, and attribute does nothing. These are really meant to be the 9 columns. Edit: There's a bug report related to this. Edit 2: I've made a pull request to clarify this and fix the aforementioned bug report. Edit 3: I realized that I ...


Only top voted, non community-wiki answers of a minimum length are eligible