Skip to content

Commit 1c66fa0

Browse files
authored
Merge pull request #67 from openvax/feature/gencode-aliases-and-versions
Add attribute_aliases + cast_version_columns to read_gtf (closes #63, #64)
2 parents d3313ac + 668dfb0 commit 1c66fa0

5 files changed

Lines changed: 716 additions & 2 deletions

File tree

gtfparse/__init__.py

Lines changed: 5 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -14,16 +14,20 @@
1414
from .create_missing_features import create_missing_features
1515
from .parsing_error import ParsingError
1616
from .read_gtf import (
17+
GENCODE_BIOTYPE_ALIASES,
18+
INTEGER_VERSION_COLUMNS,
1719
REQUIRED_COLUMNS,
1820
parse_gtf,
1921
parse_gtf_and_expand_attributes,
2022
parse_gtf_pandas,
2123
read_gtf,
2224
)
2325

24-
__version__ = "2.6.3"
26+
__version__ = "2.7.0"
2527

2628
__all__ = [
29+
"GENCODE_BIOTYPE_ALIASES",
30+
"INTEGER_VERSION_COLUMNS",
2731
"REQUIRED_COLUMNS",
2832
"ParsingError",
2933
"__version__",

gtfparse/read_gtf.py

Lines changed: 125 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -13,6 +13,7 @@
1313
import logging
1414
from os.path import exists
1515

16+
import pandas as pd
1617
import polars
1718

1819
from .attribute_parsing import expand_attribute_strings
@@ -22,6 +23,27 @@
2223
logger = logging.getLogger(__name__)
2324

2425

26+
# GENCODE GTFs use *_type where Ensembl GTFs use *_biotype. Pass this
27+
# (or a superset) as `attribute_aliases` to read_gtf to normalize a
28+
# GENCODE-format GTF onto the Ensembl column names that downstream
29+
# tools like pyensembl expect.
30+
GENCODE_BIOTYPE_ALIASES = {
31+
"gene_type": "gene_biotype",
32+
"transcript_type": "transcript_biotype",
33+
}
34+
35+
36+
# Ensembl-style attribute columns that are always integer-valued when
37+
# present. read_gtf casts these from string to pandas nullable Int64 by
38+
# default; pass cast_version_columns=False to keep them as strings.
39+
INTEGER_VERSION_COLUMNS = (
40+
"gene_version",
41+
"transcript_version",
42+
"protein_version",
43+
"exon_version",
44+
)
45+
46+
2547
"""
2648
Columns of a GTF file:
2749
@@ -183,6 +205,65 @@ def parse_gtf_and_expand_attributes(
183205
)
184206

185207

208+
def _apply_attribute_aliases(result_df, attribute_aliases):
209+
"""
210+
Rename alias attribute columns onto canonical names in-place.
211+
212+
For each (alias -> canonical) pair, in iteration order:
213+
* if only the alias is present, rename it to the canonical name.
214+
* if both are present, drop the alias and warn (canonical wins).
215+
* if neither is present, do nothing.
216+
217+
When two aliases target the same canonical (e.g. both ``gene_type``
218+
and a hypothetical ``gene_kind`` map to ``gene_biotype``), the first
219+
rename in iteration order wins; subsequent aliases targeting an
220+
already-renamed canonical are treated as collisions, dropped, and
221+
warned about.
222+
"""
223+
if not attribute_aliases:
224+
return result_df
225+
columns_present = set(result_df.columns)
226+
rename_map = {}
227+
drop_aliases = []
228+
for alias, canonical in attribute_aliases.items():
229+
if alias not in columns_present:
230+
continue
231+
if canonical in columns_present:
232+
logger.warning(
233+
"Both alias column '%s' and canonical column '%s' are present; "
234+
"dropping alias and keeping canonical values.",
235+
alias,
236+
canonical,
237+
)
238+
drop_aliases.append(alias)
239+
else:
240+
rename_map[alias] = canonical
241+
# Reflect the rename in the running column set so a later
242+
# alias mapping to the same canonical sees the collision
243+
# instead of silently producing a duplicate-named column.
244+
columns_present.discard(alias)
245+
columns_present.add(canonical)
246+
if drop_aliases:
247+
result_df = result_df.drop(columns=drop_aliases)
248+
if rename_map:
249+
result_df = result_df.rename(columns=rename_map)
250+
return result_df
251+
252+
253+
def _cast_version_columns(result_df, version_columns=INTEGER_VERSION_COLUMNS):
254+
"""
255+
Cast known Ensembl *_version attribute columns from strings to
256+
pandas nullable Int64 in-place. Missing/empty values become pd.NA.
257+
"""
258+
for column_name in version_columns:
259+
if column_name not in result_df.columns:
260+
continue
261+
result_df[column_name] = pd.to_numeric(
262+
result_df[column_name].replace("", None), errors="coerce"
263+
).astype("Int64")
264+
return result_df
265+
266+
186267
def read_gtf(
187268
filepath_or_buffer,
188269
expand_attribute_column=True,
@@ -192,6 +273,8 @@ def read_gtf(
192273
usecols=None,
193274
features=None,
194275
result_type="polars",
276+
attribute_aliases=None,
277+
cast_version_columns=True,
195278
):
196279
"""
197280
Parse a GTF into a dictionary mapping column names to sequences of values.
@@ -231,13 +314,44 @@ def read_gtf(
231314
result_type : One of 'polars', 'pandas', or 'dict'
232315
Default behavior is to return a Polars DataFrame, but will convert to
233316
Pandas DataFrame or dictionary if specified.
317+
318+
attribute_aliases : dict of str -> str, optional
319+
Maps alias attribute names onto canonical ones. After attributes
320+
are expanded into columns, each alias column is renamed to its
321+
canonical name when the canonical column is absent. If both are
322+
present the alias is dropped and a warning is logged. Pass
323+
`GENCODE_BIOTYPE_ALIASES` to normalize a GENCODE GTF's
324+
`gene_type`/`transcript_type` onto Ensembl's
325+
`gene_biotype`/`transcript_biotype`.
326+
327+
cast_version_columns : bool
328+
When True (default), cast the well-known integer version
329+
attribute columns (`gene_version`, `transcript_version`,
330+
`protein_version`, `exon_version`) from strings to pandas
331+
nullable Int64 when present. Set to False to keep them as
332+
strings.
234333
"""
235334
if type(filepath_or_buffer) is str and not exists(filepath_or_buffer):
236335
raise ValueError("GTF file does not exist: %s" % filepath_or_buffer)
237336

337+
# If usecols asks for a canonical column that's only present in the
338+
# GTF under an alias name, expand the parse-time column filter to
339+
# also pull the alias through — otherwise it gets dropped at parse
340+
# time before _apply_attribute_aliases can see it. The end-of-function
341+
# usecols filter still narrows the result down to the canonical name.
342+
parse_usecols = usecols
343+
if usecols is not None and attribute_aliases:
344+
usecols_set = set(usecols)
345+
parse_usecols = set(usecols_set)
346+
for alias, canonical in attribute_aliases.items():
347+
if canonical in usecols_set:
348+
parse_usecols.add(alias)
349+
238350
if expand_attribute_column:
239351
result_df = parse_gtf_and_expand_attributes(
240-
filepath_or_buffer, restrict_attribute_columns=usecols, features=features
352+
filepath_or_buffer,
353+
restrict_attribute_columns=parse_usecols,
354+
features=features,
241355
)
242356
else:
243357
result_df = parse_gtf(result_df, features=features)
@@ -267,6 +381,16 @@ def wrapped_fn(x):
267381
column_type = column_cast_types[column_name]
268382
result_df[column_name] = result_df[column_name].astype(column_type)
269383

384+
# Rename alias attribute columns onto their canonical names. Done before
385+
# infer_biotype_column so an aliased gene_biotype/transcript_biotype is
386+
# visible to the inference logic.
387+
result_df = _apply_attribute_aliases(result_df, attribute_aliases)
388+
389+
# Cast Ensembl *_version columns from strings to nullable integers so
390+
# downstream consumers (e.g. pyensembl) don't have to int(...) themselves.
391+
if cast_version_columns:
392+
result_df = _cast_version_columns(result_df)
393+
270394
# Hackishly infer whether the values in the 'source' column of this GTF
271395
# are actually representing a biotype by checking for the most common
272396
# gene_biotype and transcript_biotype value 'protein_coding'

tests/data/gencode.head.gtf

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,10 @@
1+
##description: GENCODE-style synthetic head for gtfparse tests
2+
##provider: gtfparse-tests
3+
##format: gtf
4+
chr1 HAVANA gene 11869 14409 . + . gene_id "ENSG00000223972"; gene_version "5"; gene_type "transcribed_unprocessed_pseudogene"; gene_name "DDX11L1";
5+
chr1 HAVANA transcript 11869 14409 . + . gene_id "ENSG00000223972"; gene_version "5"; transcript_id "ENST00000456328"; transcript_version "2"; gene_type "transcribed_unprocessed_pseudogene"; transcript_type "processed_transcript"; gene_name "DDX11L1"; transcript_name "DDX11L1-202";
6+
chr1 HAVANA exon 11869 12227 . + . gene_id "ENSG00000223972"; gene_version "5"; transcript_id "ENST00000456328"; transcript_version "2"; exon_id "ENSE00002234944"; exon_version "1"; exon_number "1"; gene_type "transcribed_unprocessed_pseudogene"; transcript_type "processed_transcript";
7+
chr1 HAVANA gene 65419 71585 . + . gene_id "ENSG00000186092"; gene_version "6"; gene_type "protein_coding"; gene_name "OR4F5";
8+
chr1 HAVANA transcript 65419 71585 . + . gene_id "ENSG00000186092"; gene_version "6"; transcript_id "ENST00000641515"; transcript_version "2"; gene_type "protein_coding"; transcript_type "protein_coding"; gene_name "OR4F5"; transcript_name "OR4F5-201";
9+
chr1 HAVANA CDS 65565 65573 . + 0 gene_id "ENSG00000186092"; gene_version "6"; transcript_id "ENST00000641515"; transcript_version "2"; exon_number "1"; protein_id "ENSP00000493376"; protein_version "2"; gene_type "protein_coding"; transcript_type "protein_coding";
10+
chr1 HAVANA exon 65565 65573 . + . gene_id "ENSG00000186092"; gene_version "6"; transcript_id "ENST00000641515"; transcript_version "2"; exon_id "ENSE00003812156"; exon_version "1"; exon_number "1"; gene_type "protein_coding"; transcript_type "protein_coding";

tests/data/gencode.real.head.gtf

Lines changed: 13 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,13 @@
1+
##description: evidence-based annotation of the human genome (GRCh38) - test excerpt
2+
##provider: GENCODE
3+
##contact: gencode-help@ebi.ac.uk
4+
##format: gtf
5+
##date: 2024-01-01
6+
chr1 HAVANA gene 11869 14409 . + . gene_id "ENSG00000223972.5"; gene_type "transcribed_unprocessed_pseudogene"; gene_name "DDX11L1"; level 2; hgnc_id "HGNC:37102"; havana_gene "OTTHUMG00000000961.2";
7+
chr1 HAVANA transcript 11869 14409 . + . gene_id "ENSG00000223972.5"; transcript_id "ENST00000456328.2"; gene_type "transcribed_unprocessed_pseudogene"; transcript_type "processed_transcript"; gene_name "DDX11L1"; transcript_name "DDX11L1-202"; level 2; transcript_support_level "1"; hgnc_id "HGNC:37102"; tag "basic"; havana_gene "OTTHUMG00000000961.2"; havana_transcript "OTTHUMT00000362751.1";
8+
chr1 HAVANA exon 11869 12227 . + . gene_id "ENSG00000223972.5"; transcript_id "ENST00000456328.2"; exon_number 1; exon_id "ENSE00002234944.1"; gene_type "transcribed_unprocessed_pseudogene"; transcript_type "processed_transcript"; gene_name "DDX11L1"; transcript_name "DDX11L1-202"; level 2; hgnc_id "HGNC:37102"; tag "basic"; havana_gene "OTTHUMG00000000961.2"; havana_transcript "OTTHUMT00000362751.1";
9+
chr1 HAVANA gene 65419 71585 . + . gene_id "ENSG00000186092.7"; gene_type "protein_coding"; gene_name "OR4F5"; level 2; hgnc_id "HGNC:14825"; havana_gene "OTTHUMG00000001094.4";
10+
chr1 HAVANA transcript 65419 71585 . + . gene_id "ENSG00000186092.7"; transcript_id "ENST00000641515.2"; gene_type "protein_coding"; transcript_type "protein_coding"; gene_name "OR4F5"; transcript_name "OR4F5-201"; level 2; tag "MANE_Select"; havana_gene "OTTHUMG00000001094.4";
11+
chr1 HAVANA CDS 65565 65573 . + 0 gene_id "ENSG00000186092.7"; transcript_id "ENST00000641515.2"; exon_number 1; protein_id "ENSP00000493376.2"; gene_type "protein_coding"; transcript_type "protein_coding"; gene_name "OR4F5"; tag "MANE_Select";
12+
chr1 HAVANA exon 65565 65573 . + . gene_id "ENSG00000186092.7"; transcript_id "ENST00000641515.2"; exon_number 1; exon_id "ENSE00003812156.1"; gene_type "protein_coding"; transcript_type "protein_coding"; gene_name "OR4F5"; tag "MANE_Select";
13+
chrM ENSEMBL gene 3307 4262 . + . gene_id "ENSG00000198888.2"; gene_type "protein_coding"; gene_name "MT-ND1"; level 3;

0 commit comments

Comments
 (0)