# How to Convert BED to GFF and GTF: A Step-by-Step Guide You can convert a BED file to GFF3 or GTF in a few minutes once you decide what each interval is supposed to represent. By the end of this guide you will have a working one-liner for BED3 and BED6 files, a Python script that rebuilds gene, mRNA, exon and CDS features from BED12, and two tested routes from GFF3 to GTF. You will also be able to validate the output and recognize the coordinate and attribute errors that break most conversions. Have these on hand before you start. A BED file with a known column count and a known meaning for each interval. The matching genome FASTA, indexed with `samtools faidx`, so you can check chromosome names. `awk` and Python 3.8 or newer. Optional but recommended: `gffread`, the AGAT toolkit, and the GenomeTools `gt` binary for validation. Install them in one environment: ```bash conda create -n bed2gff -c bioconda -c conda-forge \ gffread agat genometools-genometools samtools python=3.11 conda activate bed2gff ``` The hard part of this job is not the syntax. It is the data model. BED describes intervals and nothing else. GFF3 and GTF describe features that belong to a hierarchy, and they expect a gene ID, a transcript ID and a feature type on every line. A conversion that only shifts coordinates produces a file that parses but does not behave correctly in downstream tools such as featureCounts, HTSeq, or a genome browser track hub. ## What Each Format Encodes ### BED BED stands for Browser Extensible Data, and the format was defined by the [UCSC Genome Browser](/knowledge/bioinformatics/ucsc-genome-browser-beginner-tutorial-and-tips) group. A BED record has 3 required columns and up to 12 optional ones, tab separated: 1. `chrom`, the sequence name. 2. `chromStart`, the 0-based start coordinate. 3. `chromEnd`, the end coordinate, which is exclusive. 4. `name`, the feature label. 5. `score`, an integer from 0 to 1000. 6. `strand`, one of `+`, `-`, or `.`. 7. `thickStart`, the start of the translated or otherwise highlighted region. 8. `thickEnd`, the end of that region. 9. `itemRgb`, a display color written as `R,G,B`. 10. `blockCount`, the number of sub-blocks. 11. `blockSizes`, a comma-separated list of block lengths. 12. `blockStarts`, a comma-separated list of block offsets relative to `chromStart`. BED is 0-based and half-open. The record `chr1 100 200` covers the bases numbered 101 through 200 in 1-based counting, or the half-open interval [100, 200). The half-open convention means the length of a feature is simply `end - start`, which is convenient for arithmetic and a frequent source of confusion when moving to 1-based formats. ### GFF3 GFF3 has exactly 9 tab-separated columns: `seqid`, `source`, `type`, `start`, `end`, `score`, `strand`, `phase`, and `attributes`. It is 1-based and inclusive on both ends. The interval `start=101, end=200` in GFF3 refers to the same 100 bases as `100 200` in BED. Three things distinguish GFF3 from a flat list of intervals: - The first line must be `##gff-version 3`, and directive lines that start with `##` must appear before any feature line. - The `type` column is drawn from the Sequence Ontology, a controlled vocabulary that includes terms such as `gene`, `mRNA`, `exon`, `CDS`, `five_prime_UTR`, and `biological_region`. - The attributes column uses `tag=value` pairs, and the reserved tags `ID` and `Parent` encode the feature hierarchy. An exon points to its transcript with `Parent=tx1`, and the transcript points to its gene with `Parent=gene1`. GFF3 also supports multiple values per tag, separated by commas, and it allows URL-encoded characters for values that contain reserved characters. ### GTF GTF, sometimes written GTF2, also has 9 columns, and the coordinates are identical to GFF3: 1-based, inclusive. The difference is in the attribute column, which uses space-separated `key "value"` pairs terminated by a semicolon. Every line of a gene or transcript feature is expected to carry `gene_id` and `transcript_id`: ``` chr1 bed2gff exon 101 200 . + . gene_id "GENE1"; transcript_id "GENE1.t1"; exon_number "1"; ``` GTF has no formal hierarchy beyond that shared ID pairing. It also has no `##` directives in common use, and the `frame` column is the GTF name for what GFF3 calls `phase`. Many tools that claim to read GTF are strict about the order of the two required attributes, so keep `gene_id` before `transcript_id`. ### Coordinate systems side by side | Property | BED | GFF3 | GTF | |--|--|--|--| | Columns | 3 to 12 | 9 | 9 | | Origin | 0-based | 1-based | 1-based | | End | Exclusive | Inclusive | Inclusive | | Feature type | None (implicit) | Sequence Ontology term | Free text, conventionally SO terms | | Hierarchy | None | `ID` and `Parent` | `gene_id` and `transcript_id` | | Attribute syntax | Not applicable | `tag=value` | `key "value"` | | Header lines | `track`, `browser` (optional) | `##gff-version 3`, `##sequence-region` | Comment lines only | | Phase column | Not applicable | `0`, `1`, `2`, or `.` | `0`, `1`, `2`, or `.` | The conversion formula for the start coordinate is the same in both directions: ``` GFF_start = BED_start + 1 BED_start = GFF_start - 1 GFF_end = BED_end BED_end = GFF_end ``` Nothing else about the coordinates changes. Every other decision is a modeling decision. ## Where BED, GFF, GTF, FASTA, and FASTQ Fit in a Pipeline It helps to keep the file families straight, because the word "bioinformatics file formats" covers at least three unrelated jobs. FASTA and FASTQ hold sequence. A FASTA record is a `>` header line followed by sequence lines with no quality information. A FASTQ record is four lines: an `@` header, the sequence, a `+` separator, and a quality string of the same length encoded in Phred+33 ASCII characters. FASTQ is the raw output of an Illumina or Nanopore run, and FASTA is what you get after quality scores are no longer needed, for example when you build a reference or run a [multiple sequence alignment](/blog/guides/multiple-sequence-alignment). Annotation formats hold coordinates, not sequence. BED, GFF3 and GTF all point at positions in a FASTA file. Alignment and variant formats are a third family: SAM and BAM for reads mapped to a reference, CRAM as the compressed archival variant, and VCF or BCF for variants. The practical consequence is that a BED to GFF conversion always depends on a reference sequence whose naming matches the annotation. If your BED says `chr1` and your FASTA says `1`, the converted GFF3 will be syntactically perfect and completely unusable. ## Step-by-Step: Convert BED to GFF3 ### Step 1: Decide what each interval represents Write this down before touching the file. An interval may be a gene, a transcript, an exon, a regulatory element, a peak from a ChIP-seq experiment, or an arbitrary region of interest. The `type` column must reflect that choice, and the choice must be consistent across the file. Three common situations: - A peak file from a ChIP-seq or ATAC-seq experiment. Every interval is `biological_region` or a more specific Sequence Ontology term such as `enhancer`. There is no gene or transcript hierarchy, so each record gets a unique `ID` and no `Parent`. - A BED6 file of gene bodies. Each interval maps to `gene`. If you also have transcript models, this is where you rebuild a hierarchy. - A BED12 file exported from the UCSC Table Browser. The blocks are exons and the thick region is the CDS. This is the only BED variant that carries enough structure to build a real annotation, and it needs a script rather than a one-liner. ### Step 2: Clean and sort the BED file Real BED files arrive with `track` lines, `browser` lines, comment lines, CRLF endings, and inconsistent column counts. Strip them first. ```bash # Inspect the file head -n 3 input.bed awk -F'\t' '{print NF}' input.bed | sort -n | uniq -c # Remove track, browser, comment and blank lines, convert CRLF awk 'BEGIN{FS=OFS="\t"} /^track/ || /^browser/ || /^#/ || NF == 0 {next} NF < 3 {print "short line at " NR > "/dev/stderr"; next} {sub(/\r$/, ""); print}' input.bed > cleaned.bed # Sort by sequence then start sort -k1,1 -k2,2n cleaned.bed > cleaned.sorted.bed ``` Then compare chromosome names against the reference: ```bash samtools faidx genome.fa cut -f1 cleaned.sorted.bed | sort -u > bed.chroms cut -f1 genome.fa.fai | sort -u > fa.chroms comm -23 bed.chroms fa.chroms ``` Any line printed by that `comm` command is a chromosome that exists in your BED and not in your FASTA. Fix the naming before conversion. The usual offenders are `chr` prefixes and mitochondrial names, where UCSC uses `chrM` and Ensembl and NCBI use `MT`. ### Step 3: Convert BED3 or BED6 to GFF3 with awk This is the conversion for flat interval files with no hierarchy. ```bash { echo "##gff-version 3" awk 'BEGIN{FS=OFS="\t"} function esc(s) { gsub(/%/, "%25", s) gsub(/;/, "%3B", s) gsub(/=/, "%3D", s) gsub(/&/, "%26", s) gsub(/,/, "%2C", s) return s } { name = ($4 != "" && $4 != ".") ? $4 : "feature_" NR score = ($5 != "" && $5 != ".") ? $5 : "." strand = ($6 == "+" || $6 == "-") ? $6 : "." printf "%s\tbed2gff\tbiological_region\t%d\t%d\t%s\t%s\t.\tID=%s\n", \ $1, $2 + 1, $3, score, strand, esc(name) }' cleaned.sorted.bed } > features.gff3 ``` The `esc` function handles the five reserved characters in GFF3 attribute values. Without it, a feature name containing a semicolon truncates the attribute list and the file fails validation. Expected output for a two-line BED file: ``` ##gff-version 3 chr1 bed2gff biological_region 101 200 . + . ID=peak_1 chr1 bed2gff biological_region 501 650 12.4 - . ID=peak_2 ``` Note the start coordinate of the first feature. BED `100` became GFF3 `101`. The end stayed at `200` because BED ends are already exclusive and GFF ends are inclusive, so the numerical value is unchanged. Two refinements are worth making. First, replace the generic type with the most specific Sequence Ontology term that fits your experiment, because tools that filter by type will otherwise miss your features. Second, if the BED name is not guaranteed unique, add a prefix or a line counter to the `ID`, since duplicate IDs make a GFF3 file invalid. ### Step 4: Convert BED12 to GFF3 with a Python script BED12 carries exon blocks, a coding region, and a strand, which is enough to build a gene, mRNA, exon and CDS hierarchy. The script below does that and computes the CDS phase correctly, including on the minus strand. ```python #!/usr/bin/env python3 """Convert a BED12 file to GFF3 with gene, mRNA, exon and CDS features.""" import sys def parse_bed12(path): with open(path) as fh: for line in fh: line = line.rstrip("\n") if not line or line.startswith(("#", "track", "browser")): continue f = line.split("\t") if len(f) < 12: raise ValueError("not a BED12 line: " + line[:80]) yield { "chrom": f[0], "start": int(f[1]), "end": int(f[2]), "name": f[3], "score": f[4], "strand": f[5], "thick_start": int(f[6]), "thick_end": int(f[7]), "sizes": [int(x) for x in f[10].rstrip(",").split(",")], "starts": [int(x) for x in f[11].rstrip(",").split(",")], } def main(bed_path, out_path, source="bed2gff"): seen = set() with open(out_path, "w") as out: out.write("##gff-version 3\n") for r in parse_bed12(bed_path): name = r["name"] or "feature" if name in seen: # duplicate names break GFF3 name = name + "_" + str(r["start"]) seen.add(name) gene_id, tx_id = name, name + ".t1" strand = r["strand"] if r["strand"] in ("+", "-") else "." score = r["score"] if r["score"] not in ("", ".") else "." g_start, g_end = r["start"] + 1, r["end"] out.write(f'{r["chrom"]}\t{source}\tgene\t{g_start}\t{g_end}\t.\t{strand}\t.\tID={gene_id}\n') out.write(f'{r["chrom"]}\t{source}\tmRNA\t{g_start}\t{g_end}\t{score}\t{strand}\t.\t' f'ID={tx_id};Parent={gene_id}\n') for i, (size, rel) in enumerate(zip(r["sizes"], r["starts"]), start=1): e_start = r["start"] + rel + 1 e_end = r["start"] + rel + size out.write(f'{r["chrom"]}\t{source}\texon\t{e_start}\t{e_end}\t.\t{strand}\t.\t' f'ID={tx_id}.exon{i};Parent={tx_id}\n') # Coding blocks are the intersection of exons with the thick region cds = [] if r["thick_end"] > r["thick_start"]: for size, rel in zip(r["sizes"], r["starts"]): b_start = r["start"] + rel b_end = b_start + size c_start = max(b_start, r["thick_start"]) c_end = min(b_end, r["thick_end"]) if c_end > c_start: cds.append((c_start, c_end)) cds.sort(reverse=(strand == "-")) # translation order length_so_far = 0 for i, (c_start, c_end) in enumerate(cds, start=1): phase = (3 - (length_so_far % 3)) % 3 out.write(f'{r["chrom"]}\t{source}\tCDS\t{c_start + 1}\t{c_end}\t.\t{strand}\t{phase}\t' f'ID={tx_id}.cds{i};Parent={tx_id}\n') length_so_far += c_end - c_start if __name__ == "__main__": main(sys.argv[1], sys.argv[2]) ``` Run it and inspect the first transcript: ```bash python3 bed12_to_gff3.py cleaned.sorted.bed annotation.gff3 grep -P '\t(mRNA|exon|CDS)\t' annotation.gff3 | head -n 6 ``` A correctly converted three-exon transcript looks like this: ``` chr1 bed2gff gene 1001 1800 . + . ID=TX1 chr1 bed2gff mRNA 1001 1800 . + . ID=TX1.t1;Parent=TX1 chr1 bed2gff exon 1001 1100 . + . ID=TX1.t1.exon1;Parent=TX1.t1 chr1 bed2gff exon 1301 1400 . + . ID=TX1.t1.exon2;Parent=TX1.t1 chr1 bed2gff exon 1601 1800 . + . ID=TX1.t1.exon3;Parent=TX1.t1 chr1 bed2gff CDS 1051 1100 . + 0 ID=TX1.t1.cds1;Parent=TX1.t1 ``` The `phase` value is defined as the number of bases to remove from the start of the feature to reach the first base of the next codon. The script tracks the cumulative coding length in translation order, so the first CDS is always phase 0 and each later segment gets `(3 - (length % 3)) % 3`. For the minus strand it sorts blocks in descending coordinate order before accumulating. If your BED12 records come from the UCSC Table Browser, there is a shorter route through the UCSC utilities, which go BED12 to genePred to GTF in one pipe: ```bash bedToGenePred cleaned.sorted.bed stdout \ | genePredToGtf file stdin -source=bed2gff -utr annotation.gtf ``` The `-utr` flag writes explicit 5UTR and 3UTR features. This route is well tested and produces GTF that most count tools accept directly. It does not produce GFF3. ### Step 5: Add sequence region directives and confirm sort order Strict GFF3 consumers appreciate a `##sequence-region` line for each sequence, placed after `##gff-version 3` and before any feature. ```bash awk 'BEGIN{FS=OFS="\t"} !/^#/ {print $1, $4, $5}' annotation.gff3 \ | sort -k1,1 -k2,2n \ | awk '!seen[$1]++ {print "##sequence-region", $1, $2, $3}' > regions.txt cat regions.txt annotation.gff3 > annotation.header.gff3 ``` A GFF3 file does not have to be sorted, but most tools assume it is, and sorting by sequence then start makes diff-based checks readable. ### Step 6: Validate the GFF3 before going further The GenomeTools validator checks structure, attribute syntax, and parent references: ```bash gt gff3validator annotation.gff3 ``` A valid file produces no output and returns exit status 0. A broken parent reference, a missing `##gff-version` line, or an unescaped attribute character produces a specific error with a line number. Three quick checks catch most real problems. ```bash # Feature type counts awk -F'\t' '!/^#/ {print $3}' annotation.gff3 | sort | uniq -c # Every exon parent must exist as an mRNA ID comm -23 \ <(awk -F'\t' '!/^#/ && $3=="exon" {n=split($9,a,";"); for(i=1;i<=n;i++) if(a[i] ~ /^Parent=/){sub(/^Parent=/,"",a[i]); print a[i]}}' annotation.gff3 | sort -u) \ <(awk -F'\t' '!/^#/ && $3=="mRNA" {n=split($9,a,";"); for(i=1;i<=n;i++) if(a[i] ~ /^ID=/){sub(/^ID=/,"",a[i]); print a[i]}}' annotation.gff3 | sort -u) # Coordinate spot check: the GFF start minus one must equal the original BED start awk -F'\t' '!/^#/ && $3=="gene" {print $1, $4-1, $5}' annotation.gff3 | head head cleaned.sorted.bed ``` The two lists in the spot check should be identical once you account for the features you chose to model. ## Step-by-Step: Convert GFF3 to GTF GTF is often the final target because count tools and differential expression packages accept it more reliably than GFF3. Two converters cover almost every case. ### Route A: gffread ```bash gffread annotation.gff3 -T -o annotation.gtf head -n 3 annotation.gtf ``` The `-T` flag tells gffread to write GTF instead of GFF3. gffread translates `ID` and `Parent` into `gene_id` and `transcript_id`, and it will infer the gene ID from the transcript ID if the parent gene feature is missing. That inference is convenient and occasionally wrong, so check the output. ### Route B: AGAT ```bash agat_convert_sp_gff2gtf.pl --gff annotation.gff3 -o annotation.gtf ``` AGAT is stricter about the input hierarchy and will report missing gene or transcript parents instead of guessing. Recent AGAT releases drop the `.pl` suffix on script names, so `agat_convert_sp_gff2gtf` is the same tool. If your GFF3 has structural problems, `agat_convert_sp_gxf2gxf.pl -g annotation.gff3 -o fixed.gff3` repairs common issues such as missing parent features and misplaced directives. ### Converting back and checking the round trip A round trip is the fastest way to find attribute mapping problems. ```bash gffread annotation.gtf -o roundtrip.gff3 diff <(grep -v '^#' annotation.gff3 | cut -f1-8) \ <(grep -v '^#' roundtrip.gff3 | cut -f1-8) ``` Columns 1 through 8 should match exactly if the feature set was preserved. Differences in column 9 are normal, because attribute names and ordering are not guaranteed to survive a conversion. ### Key tools and options | Task | Command | Notes | |--|--|--| | BED3 or BED6 to GFF3 | `awk` script in Step 3 | Full control of type and ID | | BED12 to GFF3 | Python script in Step 4 | Rebuilds gene, mRNA, exon, CDS | | BED12 to GTF | `bedToGenePred` piped to `genePredToGtf` | UCSC utilities, add `-utr` for UTRs | | GFF3 to GTF | `gffread in.gff3 -T -o out.gtf` | Maps ID and Parent to gene_id and transcript_id | | GFF3 to GTF | `agat_convert_sp_gff2gtf.pl, gff in.gff3 -o out.gtf` | Strict parent checking | | GTF to GFF3 | `gffread in.gtf -o out.gff3` | Also available as `gt gtf_to_gff3` | | GFF3 to BED | `gffread in.gff3 -o out.bed` | Format chosen by output extension | | Validate GFF3 | `gt gff3validator in.gff3` | Structural and attribute checks | ## Common Mistakes and Limitations **Skipping the coordinate shift.** BED is 0-based and GFF is 1-based. A conversion that copies `chromStart` into `start` shifts every feature one base to the left. The error is invisible in a genome browser at low zoom and ruins splice site analysis at base resolution. **Treating the BED name column as a gene ID.** The `name` field is a free-text label. It may be a peak identifier, a repeat family, or a duplicate of another record. Putting it directly into `gene_id` produces GTF that runs without error and assigns reads to the wrong gene. **Leaving out `##gff-version 3`.** Many parsers reject the file, and others guess the version and guess wrong. The directive belongs on line 1, before every other line. **Mismatched chromosome naming.** `chr1` versus `1` and `chrM` versus `MT` are the most common reason a converted file fails against a reference. Check with the `comm` command in Step 2 before you convert, not after. **Unescaped attribute characters.** GFF3 requires percent-encoding of `%`, semicolon, equals, ampersand, comma, tab, and control characters inside attribute values. A gene symbol containing a comma silently splits into two values. **Wrong score semantics.** BED scores are integers from 0 to 1000, and they are often Phred-scaled or rank-based. GFF and GTF scores are floating point values or a single period. Copying a BED score of 4,312,000 into a GFF3 score column produces a file that some tools reject and others mis-scale. **Dropping BED12 blocks.** Converting BED12 with a BED3 one-liner throws away every exon boundary and the coding region. The result parses, but it is no longer the same annotation. **Duplicate IDs.** Reusing the same `ID` on two feature lines makes the GFF3 invalid. If your BED names repeat, add a coordinate or counter suffix as the Step 4 script does. **Zero-length intervals.** BED allows `start == end` for insertions and some point features. GFF3 and GTF do not represent zero-length features sensibly, and most validators flag them. Drop or pad them deliberately, and document the choice. **Assuming conversion is reversible.** The information content of BED and GFF is not the same. BED12 to GFF3 to GTF round trips lose the item color, the exact BED score, and any custom columns. GFF3 contains feature types, phase values and hierarchy that have no BED representation at all. **Forgetting that phase and frame are biological, not cosmetic.** The phase column tells the translation machinery where the next codon starts. If you set every CDS to phase 0, the annotated coding sequence will not translate correctly. The larger limitation is conceptual. BED has no vocabulary, so no converter can know whether an interval is an enhancer, an exon, or a gene. GTF has no hierarchy, so a BED12 conversion to GTF creates a flat feature set joined only by `gene_id` and `transcript_id`. Every conversion is a model of your data, and the model has to match what the downstream tool expects. If you need clinical or regulatory decisions from a specific annotation, verify the converted file against a curated release such as GENCODE or RefSeq rather than trusting a script alone. ## Frequently Asked Questions ### What is the difference between GFF and GTF? GFF3 and GTF share the same 9-column layout and the same 1-based inclusive coordinates, but GFF3 uses `tag=value` attributes with `ID` and `Parent` to encode hierarchy, while GTF uses `key "value"` pairs and relies on `gene_id` and `transcript_id` to link features. ### Why is my converted GFF start coordinate off by one? BED coordinates are 0-based and GFF3 and GTF coordinates are 1-based, so the start must be incremented by one while the end stays unchanged. Copying the BED start directly shifts every feature one base to the left. ### Can I convert BED to GFF without losing exon structure? Yes, if the source is BED12. The 12-column BED format carries block sizes and block offsets, which map directly to exon coordinates, and the thick start and thick end define the coding region that becomes CDS features. ### Does bedtools convert BED to GFF? No. bedtools reads and writes BED extensively and can read GFF as input, but it has no BED to GFF writer. Use awk, the gffread utility, the AGAT toolkit, or a short script for that step. ### How do I convert BED to GTF directly? The most reliable direct route is the UCSC pipeline: run `bedToGenePred` on a BED12 file and pipe the result into `genePredToGtf` with the `file` database argument. For BED3 and BED6 files, convert to GFF3 first and then run `gffread -T` to produce GTF. ### What is the GFF3 phase column? Phase is the number of bases to remove from the start of a CDS feature to reach the first base of the next codon, and it takes the values 0, 1, 2, or a period for non-coding features. It is computed from the cumulative coding length in translation order, which means descending coordinates on the minus strand. ### How do I validate a converted GFF3 file? Run `gt gff3validator` from GenomeTools, which checks structure, attribute syntax, and parent references. Then confirm that every exon `Parent` exists as an mRNA `ID`, and spot check that each GFF start equals the original BED start plus one. ### What is the difference between FASTQ and FASTA? FASTA stores a sequence with a header and no quality scores, while FASTQ stores the sequence plus a per-base quality string encoded in Phred+33 ASCII across four lines per record. BED, GFF3 and GTF are annotation formats that point into a FASTA reference, so they are a separate family from both.