A compilation of conversion tools for BED, SAM/BAM, psl, pslx, blast tabular and blast xml
- Post by: Joseph Hughes
- November 15, 2017
- No Comment
A wide range of formats exist for representing the comparisons of different sequences to each other: blast tabular, blast xml, psl, pslx, SAM/BAM, BED
Most of these formats can be converted from one format to another. Sometimes the format is lossless allowing for the original data to be perfectly converted
without the loss of information. Other times, the format conversion is lossy permitting the conversion of only part of the original data resulting in the loss
of some information.
Here, I have compiled the tools or UNIX commands necessary for converting from one file format to another. As you can see, I am still needing to compete some of the gaps so please let me know of any other tools which are missing.
The command is shown in full below the table.
| Conversion From Row/To Col | blast-xml | blast-tab | psl | pslx | SAM/BAM | BED |
| blast-xml | N/A | blast2tsv.xsl | blastXmlToPsl | blastXmlToPsl -pslx | blast2bam | BLAST_to_BED |
| blast-tab | perl blast2xml.pl | N/A | ?? | ?? | ?? | blast2bed |
| psl | ?? | ?? | N/A | pslToPslx | psl2sam.pl | pslToBed psl2BED |
| pslx | ?? | ?? | ?? | N/A | ?? | pslToBed |
| SAM/BAM | ?? | ?? | sam2psl.py | sam2psl.py -s | samtools view | bedtools bamtobed |
| BED | ?? | ?? | bed2psl | ?? | bedtools bedtobam | N/A |
Convert from SAM to psl using sam2psl.py
Available from: https://github.com/ndaniel/fusioncatcher/blob/master/bin/sam2psl.py
Example command:
python sam2psl.py -i test.sam -s -o test.psl
This is a lossless format conversion with the -s option , however the sequence as a read is no longer supported in the psl format.
python sam2psl.py -i test.sam -o test_no_seq.psl
Usage:
sam2psl.py [options] It takes as input a file in SAM format and it converts into a PSL format file. Options: --version show program's version number and exit -h, --help show this help message and exit -i INPUT_FILENAME, --input=INPUT_FILENAME The input file in SAM format. -4, --skip-conversion-cigar-1.3 By default if the CIGAR strings in the input SAM file are in the format defined in SAM version 1.4 (i.e. there are 'X' and '=') then the CIGAR string will be first converted into CIGAR string, which is described in SAM version 1.3, (i.e. there are no 'X' and '=' which are replaced with 'M') and afterwards into PSL format. Default is 'False'. -s, --read-seq It adds to the PSL output as column 22, the sequence of the read. This is not anymore a valid PSL format. -r REPLACE_READS_IDS, --replace-read-ids=REPLACE_READS_IDS In the reads ids (also known as query name in PSL) the string specified here will be replaced with '/' (which is used in Solexa for /1 and /2). -o OUTPUT_FILENAME, --output=OUTPUT_FILENAME The output file in PSL format.
Convert from psl to SAM
Available from the samtools legacy scripts: https://github.com/lh3/samtools-legacy/blob/master/misc/psl2sam.pl
Example command:
psl2sam.pl test.psl
This ends up being a lossy conversion as the read sequence is not in the output.
Usage
psl2sam.pl Usage: psl2sam.pl [-a 1] [-b 3] [-q 5] [-r 2]
The options are used to calculate a blast like scoring see post:
https://www.biostars.org/p/44281/
Convert psl to pslx
Using https://github.com/ENCODE-DCC/kentUtils/tree/master/bin/linux.x86_64
pslToPslx test_no_seq.psl test.fa ref.fa test.pslx
This is a lossless conversion. For usage:
pslToPslx - Convert from psl to pslx format, which includes sequences usage: pslToPslx [options] in.psl qSeqSpec tSeqSpec out.pslx qSeqSpec and tSeqSpec can be nib directory, a 2bit file, or a FASTA file. FASTA files should end in .fa, .fa.gz, .fa.Z, or .fa.bz2 and are read into memory. Options: -masked - if specified, repeats are in lower case cases, otherwise entire sequence is loader case.
Convert SAM to fasta
awk '$1~!/^@/ {print ">"$1"\n"$10}' test.sam > test.fa
Convert psl to BED
Option 1:
Using pslToBed from https://github.com/ENCODE-DCC/kentUtils/tree/master/bin/linux.x86_64
This is a lossless conversion as the standard psl doesn’t have the sequence and so the bed file doesn’t either.
pslToBed test_no_seq.psl test.bed
Usage:
bedToPsl - convert bed format files to psl format
usage:
bedToPsl chromSizes bedFile pslFile
Convert a BED file to a PSL file. This the result is an alignment.
It is intended to allow processing by tools that operate on PSL.
If the BED has at least 12 columns, then a PSL with blocks is created.
Otherwise single-exon PSLs are created.
Options:
-keepQuery - instead of creating a fake query, create PSL with identical query and
target specs. Useful if bed features are to be lifted with pslMap and one
wants to keep the source location in the lift result.
(khmerEnv)hugh01j@Alpha:~/test_format_convert$ kentUtils/bin/linux.x86_64/pslToBed
pslToBed: tranform a psl format file to a bed format file.
usage:
pslToBed psl bed
options:
-cds=cdsFile
cdsFile specifies a input cds tab-separated file which contains
genbank-style CDS records showing cdsStart..cdsEnd
e.g. NM_123456 34..305
These coordinates are assumed to be in the query coordinate system
of the psl, like those that are created from genePredToFakePsl
-posName
changes the qName field to qName:qStart-qEnd
(can be used to create links to query position on details page)
Option 2 as suggested by Alex Reynolds:
Using psl2bed from http://bedops.readthedocs.io/en/latest/content/reference/file-management/conversion/psl2bed.html
This is also lossless when used with –keep-header:
Example:
psl2bed < in.psl > out.bed
As a bonus, it uses sort-bed to make a sorted BED file, so that it is ready to use with bedops, bedmap, etc.
Usage:
convert2bed -i psl
version: 2.4.29 (typical)
author: Alex Reynolds
Converts 0-based, half-open [a-1, b) headered or headerless PSL
input into 0-based, half-open [a-1, b) extended BED or BEDOPS Starch
$ psl2bed < foo.psl > sorted-foo.psl.bed
$ psl2starch < foo.psl > sorted-foo.psl.starch
Or:
$ convert2bed -i psl < foo.psl > sorted-foo.psl.bed
$ convert2bed -i psl -o starch < foo.psl > sorted-foo.psl.starch
We make no assumptions about sort order from converted output. Apply
the usage case displayed to pass data to the BEDOPS sort-bed application
which generates lexicographically-sorted BED data as output.
If you want to skip sorting, use the --do-not-sort option:
$ psl2bed --do-not-sort < foo.psl > unsorted-foo.psl.bed
The PSL specification (http://genome.ucsc.edu/goldenPath/help/blatSpec.html)
contains 21 columns, some which map to UCSC BED columns and some which do not.
PSL input can contain a header or be headerless, if the BLAT search was
performed with the -noHead option. This program can accept input in
either format.
If input is headered, you can use the --keep-header option to preserve the header
data as pseudo-BED elements that use the "_header" chromosome name. We expect this
should not cause any collision problems since PSL data should use UCSC chromosome
naming conventions.
We describe below how we map columns to BED, so that BLAT results can be losslessly
transformed back into PSL format with a simple awk statement or other similar
command that permutes columns into PSL-ordering.
We map the following PSL columns to their equivalent BED column, as follows:
- tName <--> chromosome
- tStart <--> start
- tEnd <--> stop
- qName <--> id
- qSize <--> score
- strand <--> strand
Remaining PSL columns are mapped, in order, to columns 7 through 21 in the
BED output:
- matches
- misMatches
- repMatches
- nCount
- qNumInsert
- qBaseInsert
- tNumInsert
- tBaseInsert
- qStart
- qEnd
- tSize
- blockCount
- blockSizes
- qStarts
- tStarts
PSL conversion options:
--keep-header (-k)
Preserve header section as pseudo-BED elements (requires --headered)
--split (-s)
Split record into multiple BED elements, based on tStarts field value
Other processing options:
--do-not-sort (-d)
Do not sort BED output with sort-bed (not compatible with --output=starch)
--max-mem=<value> (-m <val>)
Sets aside <value> memory for sorting BED output. For example, <value> can
be 8G, 8000M or 8000000000 to specify 8 GB of memory (default is 2G)
--sort-tmpdir=<dir> (-r <dir>)
Optionally sets [dir] as temporary directory for sort data, when used in
conjunction with --max-mem=[value], instead of the host's operating system
default temporary directory
--starch-bzip2 (-z)
Used with --output=starch, the compressed output explicitly applies the bzip2
algorithm to compress intermediate data (default is bzip2)
--starch-gzip (-g)
Used with --output=starch, the compressed output applies gzip compression on
intermediate data
--starch-note="xyz..." (-e "xyz...")
Used with --output=starch, this adds a note to the Starch archive metadata
--help | --help[-bam|-gff|-gtf|-gvf|-psl|-rmsk|-sam|-vcf|-wig] (-h | -h <fmt>)
Show general help message (or detailed help for a specified input format)
--version (-w)
Show application version
Convert pslx to BED
Using https://github.com/ENCODE-DCC/kentUtils/tree/master/bin/linux.x86_64
This is a lossy conversion as the sequence is lost
pslToBed test.pslx test.bed
Usage:
pslToBed: tranform a psl format file to a bed format file.
usage:
pslToBed psl bed
options:
-cds=cdsFile
cdsFile specifies a input cds tab-separated file which contains
genbank-style CDS records showing cdsStart..cdsEnd
e.g. NM_123456 34..305
These coordinates are assumed to be in the query coordinate system
of the psl, like those that are created from genePredToFakePsl
-posName
changes the qName field to qName:qStart-qEnd
(can be used to create links to query position on details page)
Convert BAM to BED
Using bedtools: http://bedtools.readthedocs.io/en/latest/content/tools/bamtobed.html
This is a lossy conversion as the sequence data is lost.
bedtools bamtobed -i test.bam > test_bamtobed.bed
Usage:
Tool: bedtools bamtobed (aka bamToBed) Version: v2.22.1 Summary: Converts BAM alignments to BED6 or BEDPE format. Usage: bedtools bamtobed [OPTIONS] -i <bam> Options: -bedpe Write BEDPE format. - Requires BAM to be grouped or sorted by query. -mate1 When writing BEDPE (-bedpe) format, always report mate one as the first BEDPE "block". -bed12 Write "blocked" BED format (aka "BED12"). Forces -split. http://genome-test.cse.ucsc.edu/FAQ/FAQformat#format1 -split Report "split" BAM alignments as separate BED entries. Splits only on N CIGAR operations. -splitD Split alignments based on N and D CIGAR operators. Forces -split. -ed Use BAM edit distance (NM tag) for BED score. - Default for BED is to use mapping quality. - Default for BEDPE is to use the minimum of the two mapping qualities for the pair. - When -ed is used with -bedpe, the total edit distance from the two mates is reported. -tag Use other NUMERIC BAM alignment tag for BED score. - Default for BED is to use mapping quality. Disallowed with BEDPE output. -color An R,G,B string for the color used with BED12 format. Default is (255,0,0). -cigar Add the CIGAR string to the BED entry as a 7th column.
Convert BED to BAM
Create the genome file for bed
samtools faidx ref.fa
awk -v OFS='\t' {'print $1,$2'} ref.fai > ref.txt
Using the genome file and BED file to produce the BAM file.
bedtools bedtobam -i test_bamtobed.bed -g ref_revcomp.txt > test_bedtobam.bam
Usage:
Tool: bedtools bedtobam (aka bedToBam) Version: v2.22.1 Summary: Converts feature records to BAM format. Usage: bedtools bedtobam [OPTIONS] -i <bed/gff/vcf> -g <genome> Options: -mapq Set the mappinq quality for the BAM records. (INT) Default: 255 -bed12 The BED file is in BED12 format. The BAM CIGAR string will reflect BED "blocks". -ubam Write uncompressed BAM output. Default writes compressed BAM. Notes: (1) BED files must be at least BED4 to create BAM (needs name field).
The sequence is not present in the BED file so is absent from BAM as well. This is a lossy format conversion.
Additionally, there are differences in the number of read compared to the original file.
Convert BED to psl
Using https://github.com/ENCODE-DCC/kentUtils/tree/master/bin/linux.x86_64
This is a lossless conversion as neither the BED nor psl contain sequence information
bedToPsl Longest_revcomp.txt test.bed test_bedtopsl.psl
Usage:
bedToPsl - convert bed format files to psl format usage: bedToPsl chromSizes bedFile pslFile Convert a BED file to a PSL file. This the result is an alignment. It is intended to allow processing by tools that operate on PSL. If the BED has at least 12 columns, then a PSL with blocks is created. Otherwise single-exon PSLs are created. Options: -keepQuery - instead of creating a fake query, create PSL with identical query and target specs. Useful if bed features are to be lifted with pslMap and one wants to keep the source location in the lift result.
Preparing blast-xml format
makeblastdb -dbtype nucl -in Longest_revcomp.fa blastn -query test.fa -db Longest_revcomp.fa -out test.blastxml -outfmt 5
Preparing blast-tab format
blastn -query test.fa -db Longest_revcomp.fa -out test.blasttab -outfmt 6
Blast-xml to psl
Using https://github.com/ENCODE-DCC/kentUtils/tree/master/bin/linux.x86_64
This is a lossy conversion
blastXmlToPsl test.blastxml test_blastxmltopsl.psl
However, if you use the -pslx option, you can get lossless conversion
blastXmlToPsl -pslx test.blastxml test_blastxmltopsl.pslx
Usage:
blastXmlToPsl - convert blast XML output to PSLs usage: blastXmlToPsl [options] blastXml psl options: -scores=file - Write score information to this file. Format is: strands qName qStart qEnd tName tStart tEnd bitscore eVal qDef tDef -verbose=n - n >= 3 prints each line of file after parsing. n >= 4 dumps the result of each query -eVal=n n is e-value threshold to filter results. Format can be either an integer, double or 1e-10. Default is no filter. -pslx - create PSLX output (includes sequences for blocks) -convertToNucCoords - convert protein to nucleic alignments to nucleic to nucleic coordinates -qName=src - define element used to obtain the qName. The following values are support: o query-ID - use contents of the element if it exists, otherwise use o query-def0 - use the first white-space separated word of the element if it exists, otherwise the first word of . Default is query-def0. -tName=src - define element used to obtain the tName. The following values are support: o Hit_id - use contents of the element. o Hit_def0 - use the first white-space separated word of the element. o Hit_accession - contents of the element. Default is Hit-def0. -forcePsiBlast - treat as output of PSI-BLAST. blast-2.2.16 and maybe others indentify psiblast as blastp. Output only results of last round from PSI BLAST
Converting from blast-xml to SAM/BAM
Using https://github.com/guyduche/Blast2Bam
This is a lossless conversion with sequence and read quality introduced.
blast2bam -o test_blastxmltosam.sam test.blastxml ref.fa reads_1.fq reads_2.fq
Usage
Blast2Bam. Last compilation: Jun 27 2017 at 15:21:50.
Usage: blast2bam [options] [FastQ_2]
Options:
--output | -o FILE Output file (default: stdout)
--interleaved | -p Interleaved data
--readGroup | -R STR Read group header line '@RG\tID:foo'
--minAlignLength | -W INT Discard alignments shorter than [INT]
--shortCigar | -c Short version of the CIGAR string ('M' instead of '=' and 'X')
--posOnChr | -z Adjust the alignment position to the first position of the reference
--help | -h Get help (this screen)
Subsequently converted to BAM using samtools
samtools view -b test_blastxmltosam.sam > test_blastxmltosam.bam
Blast-xml to BED
Using https://github.com/mojaveazure/BLAST_to_BED
Command
BLAST_to_BED.py -x test.blastxml -o test_blastxmltobed.bed
This is a lossy conversion as the sequence information is lost.
usage: BLAST_to_BED.py (-f FASTA FILE | -x XML FILE)
[-b BED FILE | -o OUTPUT FILE]
[-d REFERENCE BLAST DATABASE] [-e E-VALUE THRESHOLD]
[-s MAX HITS] [-m MAX HSPS] [-k]
[-X [EXCLUDE CHROMOSOMES [EXCLUDE CHROMOSOMES ...]] |
-K [KEEP ONLY CHROMOSOMES [KEEP ONLY CHROMOSOMES ...]]]
Input options:
Specify the type of input we're working with, '-f | --fasta' requires BLASTn from NCBI to be installed
-f FASTA FILE, --fasta FASTA FILE
Input FASTA file to run BLAST on, incompatible with
'-x | --xml'
-x XML FILE, --xml XML FILE
Input BLAST XML file to turn into BED file,
incompatible with '-f | --fasta'
Output options:
Specify how we're writing our output files, create new file or append to preexisting BED file
-b BED FILE, --bed BED FILE
BED file to append results to, incompatible with '-o |
--outfile'
-o OUTPUT FILE, --outfile OUTPUT FILE
Name of output file to write to, defaults to
'/home1/hugh01j/test_format_convert/output.bed',
incompatible with '-b | --bed'
BLAST options:
Options only used when BLASTING, no not need to provide with '-x | --xml'
-d REFERENCE BLAST DATABASE, --database REFERENCE BLAST DATABASE
Reference BLAST database in nucleotide format, used
only when running BLAST
-e E-VALUE THRESHOLD, --evalue E-VALUE THRESHOLD
Evalue threshold for BLAST, defaults to '1e-1', used
only when running BLAST
-s MAX HITS, --max-hits MAX HITS
Maximum hits per query, defaults to '1', used only
when running BLAST
-m MAX HSPS, --max-hsps MAX HSPS
Maximum HSPs per hit, defaults to '1', used only when
running BLAST
-k, --keep-xml Do we keep the XML results? pas '-k | --keep-xml' to
say 'yes', used only when running BLAST
Filtering options:
Options to filter the resulting BED file, can specify as many or as few chromosome/contig names as you want. You may choose to either exclude or keep, not both. Note: this only affects the immediate output of this script, we will not filter a preexisting BED file
-X [EXCLUDE CHROMOSOMES [EXCLUDE CHROMOSOMES ...]], --exclude-chrom [EXCLUDE CHROMOSOMES [EXCLUDE CHROMOSOMES ...]]
Chromosomes/Contigs to exclude from BED file,
everything else is kept, incompatible with '-K |
--keep-chrom'
-K [KEEP ONLY CHROMOSOMES [KEEP ONLY CHROMOSOMES ...]], --keep-chrom [KEEP ONLY CHROMOSOMES [KEEP ONLY CHROMOSOMES ...]]
Chromosomes/Contigs to keep from BED file, everything
else is excluded, incompatible with '-X | --exclude-
chrom'
Converting Blast tabular to blast-xml
Not completely possible due to missing information such as the alignment but see post https://www.biostars.org/p/7981/#7983
Or using the script blast2xml.pl from the blast2go google group: https://11302289526311073549.googlegroups.com/attach/ed2c446e1b1852a9/blast2xml.pl?part=0.1&view=1&vt=ANaJVrEJYYa7SZC-uvOtoKb6932qlMJWltc2p_5GrTK5Wi7jo-hw14zFroKEcLhdNcJUcQweoUJOuXk2H7wQB5q6mzDTTn211hC2OvwiWw0b5PZev-HQ7Qg
Command:
perl blast2xml.pl -i test.blastxml -o test_blasttoblastxml
Usage
perl blast2xml.pl - i|input : path of the input file (must be text blast file output) - o|output : path of the output file (by default, the same as the input file) - s|sequences : number of sequences by xml file (default inf) - hit : number of hit to print for each sequences (default inf) - hsp : number of hsp to print for each hit (default inf) - help|h|? : print this help and exit
Convert blast-xml to blast tabular
Several approaches have been suggested here https://www.biostars.org/p/7290/ but the most straight forward I have found is using the style sheet blast2tsv.xml from here:
https://github.com/lindenb/xslt-sandbox/blob/master/stylesheets/bio/ncbi/blast2tsv.xsl. This is a lossless conversion but is not a standardly formatted blast-tabular output as it contains the sequence and the aligned site information in the last two columns.
Command:
xsltproc --novalid blast2tsv.xsl test.blastxml
Blast tabular to BED
Using https://github.com/nterhoeven/blast2bed
blast2bed test.blasttab
The output will be in test.blasttab.bed, this is a lossless conversion neither blast-tabular nor BED have the sequence.
Usage
Usage: ./blast2bed The blast file should be in blast outfmt 6 or 7. See Readme.org for more details.
Converting SAM to BAM
Using samtools http://quinlanlab.org/tutorials/samtools/samtools.html
This is a lossless format conversion.
Command
samtools view -bS test.sam > test.bam
Converting BAM to SAM using samtools
This is a lossless format conversion.
samtools view -h -o test.sam test.bam
PS: I dedicate this tutorial to Sej, a great bioinformatician and friend 😉