Version Compatibility

Reference examples tested with: BioPython 1.83+, GEOparse 2.0+, R Bioconductor GEOquery 2.70+, pandas 2.2+

Before using code patterns, verify installed versions match. If versions differ:

• Python: pip show biopython geoparse then introspect signatures
• R: packageVersion('GEOquery')

If the GSE structure doesn't match expectations (missing fields, malformed series matrix), re-fetch from FTP directly and inspect the SOFT or MINiML file as source of truth.

GEO Data

"Pull expression data from GEO accession GSE..." -> GEO stores Series (GSE), Samples (GSM), Platforms (GPL), and curated DataSets (GDS, frozen 2018). The single most consequential decision is processed (series matrix) vs raw (supplementary files / linked SRA) — the answer turns on how much trust the submitter's normalization deserves.

The single most-missed gotcha: SuperSeries. A GSE may be a meta-container (!Series_relation = SuperSeries of: GSExxxxx) holding multiple sub-studies on different platforms. Naively pulling samples from a SuperSeries gives mixed Affymetrix + Illumina + RNA-seq, mis-batched.

• Python: Entrez.esearch(db='gds'), GEOparse for full series download
• R: GEOquery::getGEO() (Bioconductor; more mature than GEOparse)
• CLI: wget from ftp.ncbi.nlm.nih.gov/geo/series/...

Required Setup

pip install biopython GEOparse pandas
# OR for R-side:
# R: BiocManager::install('GEOquery')

from Bio import Entrez
Entrez.email = '[email protected]'
Entrez.api_key = 'optional'

GEO record taxonomy

GDS is dead-as-format: NCBI stopped creating new GDS records in 2018. Existing GDS still queryable but use GSE for anything current.

The SuperSeries trap

A SuperSeries (GSE) wraps multiple SubSeries, often with different platforms. Detection:

# Read the !Series_relation field from SOFT format
from Bio import Entrez
h = Entrez.esummary(db='gds', id='200122288')   # example
r = Entrez.read(h)[0]; h.close()
print(r.get('summary'))   # may or may not flag SuperSeries
# Definitive check: download SOFT and grep:
#   curl ftp://ftp.ncbi.nlm.nih.gov/geo/series/GSE122nnn/GSE122288/soft/GSE122288_family.soft.gz | zgrep Series_relation

A SuperSeries of: GSE12345 line means the SuperSeries' samples are the union of all SubSeries — almost certainly mixed-platform / mixed-batch. Process each SubSeries independently.

Symmetric trap: a paper may cite a SubSeries (SubSeries of: GSEsuper) where the wider context is essential — check both directions.

Decision matrix: processed vs raw vs SRA

Default to raw whenever possible. For Affymetrix: CEL + locally-run RMA is far more reliable than the submitter's "normalized" matrix. For RNA-seq: SRA FASTQ + locally-run alignment/quantification is the only reproducible path; submitter counts often use a private pipeline.

Series matrix files

A series matrix (GSE12345seriesmatrix.txt.gz) is a header (sample metadata as !Sample_* lines) plus a sample-by-feature expression table. The format is fragile and the values' provenance is whatever the submitter chose. Critical caveats:

• For Affymetrix: the matrix is usually RMA-normalized but submitters sometimes apply additional transforms (log2, scaling, batch correction).
• For RNA-seq: the matrix is sometimes log-CPM, sometimes raw counts, sometimes VST/rlog — read !Seriesoveralldesign and !Sampledataprocessing to know.
• The header has !Samplecharacteristicsch1 rows that hold the metadata of interest — these are submitter-formatted strings, often inconsistent within one series.

SOFT vs MINiML

Both contain the same content. SOFT is the legacy, MINiML the XML successor. GEOparse handles SOFT well; for very large series (1000+ samples) MINiML's XML structure is slower to parse.

GEOparse vs GEOquery

For production GEO workflows in R, GEOquery is the stable choice. For Python, GEOparse is the only option but verify file counts after download.

GEOmetadb status

GEOmetadb (Zhu 2008) was a SQLite mirror of GEO metadata enabling fast SQL queries. Unmaintained since 2020; downloads still work but data is stale. Modern replacement: pysradb (pysradb gsetosrp, pysradb metadata) covers most of the GEO->SRA mapping; for full GEO queries fall back to Entrez gds.

ArrayExpress -> BioStudies migration (2020)

ArrayExpress (EMBL-EBI's microarray archive, mirroring GEO) was migrated into BioStudies in 2020. Old E-MTAB-#### accessions still resolve but the API moved:

For new workflows, use BioStudies. For legacy ArrayExpress URLs in old papers, redirect via BioStudies.

Code patterns

Search GEO for studies matching a query

Goal: Find GSE accessions matching keywords + organism + study type.

Approach: ESearch on gds db with field-qualified terms; filter to gse[Entry Type]; summarize with ESummary.

Reference (BioPython 1.83+):

from Bio import Entrez
import time

Entrez.email = '[email protected]'


def search_geo(term, study_type='gse', organism=None, max_results=50):
    full_term = f'{term} AND {study_type}[Entry Type]'
    if organism:
        full_term += f' AND {organism}[Organism]'
    h = Entrez.esearch(db='gds', term=full_term, retmax=max_results)
    s = Entrez.read(h); h.close()
    if not s['IdList']:
        return []
    h = Entrez.esummary(db='gds', id=','.join(s['IdList']))
    summaries = Entrez.read(h); h.close()
    return summaries


for s in search_geo('breast cancer RNA-seq', organism='Homo sapiens', max_results=10):
    # Surface SuperSeries
    relation = s.get('summary', '')
    is_super = 'SuperSeries' in str(relation)
    print(f"  {s['Accession']:12} {s['n_samples']:>4} samples  {'[SuperSeries]' if is_super else '':12}  {s['title'][:60]}")

Detect SuperSeries before pulling data

Goal: Avoid mixing platforms by detecting SuperSeries structure first.

Approach: Download SOFT family file and read !Series_relation keys.

import gzip
import urllib.request


def check_super_or_sub_series(gse):
    prefix = gse[:-3] + 'nnn'
    url = f'https://ftp.ncbi.nlm.nih.gov/geo/series/{prefix}/{gse}/soft/{gse}_family.soft.gz'
    urllib.request.urlretrieve(url, f'{gse}.soft.gz')
    super_of = []
    sub_of = None
    with gzip.open(f'{gse}.soft.gz', 'rt') as f:
        for line in f:
            if line.startswith('!Series_relation'):
                if 'SuperSeries of' in line:
                    super_of.append(line.split('SuperSeries of: ')[1].strip())
                elif 'SubSeries of' in line:
                    sub_of = line.split('SubSeries of: ')[1].strip()
            if line.startswith('^SAMPLE'):
                break   # Speed: don't read past header
    return {'super_of': super_of, 'sub_of': sub_of}


print(check_super_or_sub_series('GSE122288'))
# {'super_of': ['GSExxxxx', 'GSEyyyyy'], 'sub_of': None}  -> SuperSeries; process subseries separately

Download series matrix with submitter caveat

import gzip
import pandas as pd


def download_series_matrix(gse):
    prefix = gse[:-3] + 'nnn'
    url = f'https://ftp.ncbi.nlm.nih.gov/geo/series/{prefix}/{gse}/matrix/{gse}_series_matrix.txt.gz'
    urllib.request.urlretrieve(url, f'{gse}_matrix.txt.gz')
    return f'{gse}_matrix.txt.gz'


def parse_series_matrix(path):
    metadata = {}
    with gzip.open(path, 'rt') as f:
        for line in f:
            if line.startswith('!series_matrix_table_begin'):
                break
            if line.startswith('!'):
                key, *vals = line.rstrip('\n').split('\t')
                metadata[key] = [v.strip('"') for v in vals]
        expr = pd.read_csv(f, sep='\t', index_col=0, comment='!')
    # Series matrix values are whatever submitter chose -- check metadata['!Sample_data_processing']
    return metadata, expr


meta, expr = parse_series_matrix(download_series_matrix('GSE123456'))
print('Sample-level data processing notes:')
for note in set(meta.get('!Sample_data_processing', [])):
    print(f'  - {note}')

Link GEO Series to SRA runs (preferred path: pysradb)

from pysradb import SRAweb


def gse_to_srr(gse):
    db = SRAweb()
    srp_df = db.gse_to_srp(gse)
    if srp_df.empty:
        return []
    srp = srp_df['study_accession'].iloc[0]
    srr_df = db.srp_to_srr(srp)
    return srr_df['run_accession'].tolist()


srrs = gse_to_srr('GSE123456')
print(f'GSE123456 -> {len(srrs)} SRR runs')

GEOparse: full Series download

import GEOparse


def get_gse(gse_id, dest='./geo_cache'):
    gse = GEOparse.get_GEO(geo=gse_id, destdir=dest)
    print(f'{gse_id}: {len(gse.gsms)} samples, {len(gse.gpls)} platforms')
    for gsm_name, gsm in list(gse.gsms.items())[:3]:
        print(f'  {gsm_name}: {gsm.metadata.get("title", ["?"])[0]}')
    return gse


# Supplementary files (raw data) -- verify file count manually after
gse = get_gse('GSE123456')
gse.download_supplementary_files(directory='./geo_cache')

R: GEOquery (more reliable supplementary download)

# Reference: Bioconductor GEOquery 2.70+ | Verify API if version differs
library(GEOquery)

gse <- getGEO('GSE123456', GSEMatrix = TRUE)
length(gse)             # one ExpressionSet per platform
head(pData(gse[[1]]))   # sample metadata
head(exprs(gse[[1]]))   # expression matrix (submitter-normalized -- verify processing notes)

# Raw / supplementary files
supp_dir <- getGEOSuppFiles('GSE123456', baseDir = './geo_cache')
list.files(rownames(supp_dir))

Find datasets by PubMed citation

def geo_from_pubmed(pmid):
    h = Entrez.elink(dbfrom='pubmed', db='gds', id=pmid)
    r = Entrez.read(h); h.close()
    if not r[0]['LinkSetDb']:
        return []
    gds_ids = [l['Id'] for l in r[0]['LinkSetDb'][0]['Link']]
    h = Entrez.esummary(db='gds', id=','.join(gds_ids))
    summaries = Entrez.read(h); h.close()
    return summaries

Failure modes

SuperSeries pulled as one experiment

• Trigger: GSE accession from a paper; turns out to be a SuperSeries wrapping multiple platforms.
• Mechanism: Default download merges all samples without flagging the structure.
• Symptom: Downstream batch correction can't recover the mixed-platform structure; spurious "batch" effects.
• Fix: Always check !Series_relation in SOFT before pulling; process SubSeries independently.

Series matrix is not what it appears to be

• Trigger: Series matrix downloaded; treated as RMA-normalized when submitter applied additional transforms.
• Mechanism: Series matrix contents are at submitter's discretion.
• Symptom: Re-analysis gives different answers than the published paper.
• Fix: Read !Sampledataprocessing to know what's in the matrix; re-normalize from raw if in doubt.

Submitter-provided RNA-seq counts mis-trusted

• Trigger: Using a *_counts.txt.gz supplementary file as the count matrix.
• Mechanism: Submitter's pipeline (aligner, GTF version, counting strategy) is rarely documented.
• Symptom: Counts don't agree with re-quantification from SRA FASTQ.
• Fix: Pull SRA FASTQ + re-quantify with a known pipeline (Salmon, kallisto, STAR + featureCounts).

Platform GPL mismatch

• Trigger: One GSE with multiple platforms; series matrix split across multiple files.
• Mechanism: GSEseriesmatrix.txt.gz is the merged one; per-platform are GSE-GPLxxxseriesmatrix.txt.gz.
• Symptom: "Missing samples" or NaN-heavy expression matrix.
• Fix: Download per-platform matrix files; check !Seriesplatformid count.

GEOparse supplementary files flakey

• Trigger: gse.downloadsupplementaryfiles() silently misses files.
• Mechanism: Known issue with the GEOparse FTP enumeration since ~2022.
• Symptom: Local cache missing CEL or counts files.
• Fix: Use R GEOquery or direct FTP wget -r on the suppl/ subdirectory.

ArrayExpress URL rot

• Trigger: Old paper links https://www.ebi.ac.uk/arrayexpress/experiments/E-MTAB-1234/.
• Mechanism: ArrayExpress migrated to BioStudies in 2020.
• Symptom: 404 or redirect.
• Fix: Use https://www.ebi.ac.uk/biostudies/arrayexpress/studies/E-MTAB-1234.

GEOmetadb stale

• Trigger: Old pipeline downloads GEOmetadb.sqlite for fast queries.
• Mechanism: GEOmetadb unmaintained since 2020.
• Symptom: Missing recent series; outdated annotations.
• Fix: Switch to pysradb for SRA-linked queries; Entrez gds for full GEO.

Common errors

References

• Edgar R, Domrachev M, Lash AE. (2002) Gene Expression Omnibus: NCBI gene expression and hybridization array data repository. Nucleic Acids Res 30:207-210.
• Barrett T, Wilhite SE, Ledoux P, et al. (2013) NCBI GEO: archive for functional genomics data sets - update. Nucleic Acids Res 41:D991-D995.
• Davis S, Meltzer PS. (2007) GEOquery: a bridge between the Gene Expression Omnibus (GEO) and BioConductor. Bioinformatics 23:1846-1847.
• Gumienny R. GEOparse: Python library to parse GEO databases. https://github.com/guma44/GEOparse (no journal publication).
• Sarkans U, Gostev M, Athar A, et al. (2018) The BioStudies database--one stop shop for all data supporting a life sciences study. Nucleic Acids Res 46:D1266-D1270.
• Lachmann A, Torre D, Keenan AB, et al. (2018) Massive mining of publicly available RNA-seq data from human and mouse. Nat Commun 9:1366. (ARCHS4)
• Wilks C, Zheng SC, Chen FY, et al. (2021) recount3: summaries and queries for large-scale RNA-seq expression and splicing. Genome Biol 22:323.

Related Skills

• entrez-search - General gds search
• entrez-link - gds <-> pubmed, bioproject links (gds->sra ELink is unreliable; use pysradb)
• sra-data - Download raw FASTQ from GEO-linked SRA runs
• expression-matrix/normalization - Re-normalize raw expression data
• rna-quantification/alignment-free-quant - Salmon/kallisto re-quantification of GEO/SRA data
• ensembl-rest - Cross-reference Ensembl IDs in series-matrix files