Method

How genEra dates a gene

genEra (Barrera-Redondo et al. 2023) assigns every protein in a genome to a phylostratum — the oldest node on the species' own lineage at which the gene can still be detected. The idea is simple enough to describe in one paragraph. The interesting part is where it breaks, and it breaks in ways that matter specifically for de novo gene claims.

The pipeline

Everything below is one focal species at a time. The ladder is the focal species' own NCBI lineage, from cellular organisms at the root down to the species itself — 19 or 20 rungs for a plant.

focal proteome 39,522 proteins · taxid DIAMOND vs NR + custom panel · ~10 h hit → taxid → lineage NCBI taxonomy project onto focal ladder deepest shared rung per hit focal lineage · oldest at top hits landing here PS 1 · cellular organisms PS 2 · Eukaryota PS 3 · Viridiplantae PS 6 · Embryophyta PS 9 · Spermatophyta PS 12 · Pentapetalae PS 16 · family PS 17 · tribe PS 19 · focal lineage 3 · moss, fern 17 · rice, pine… 64 · eudicots 88 · willows 2 · self most distant hit sets the age gene age = PS 6, Embryophyta taxonomic representativeness 100% of ranks 6→19 occupied · assign founder events (families, not genes) 4,191 @ PS 1 428 @ PS 19 · from 450 genes
Why step 5 is where the trouble lives

The age is a maximum over hits. A maximum has no tolerance for a single bad observation — and NR contains bad observations. genEra knows this, which is why step 6 exists. The next section is a real case where step 6 did its job.

The representativeness gate, on a real gene

COR15A is an Arabidopsis cold-regulated LEA protein, genuinely restricted to Brassicaceae. In this project's genEra run it came back not as "Brassicaceae" but as possible contamination or HGT, with a representativeness of 22%. Here is why that was the right answer.

DIAMOND hits for COR15A 96 hits · Brassicaceae · up to 215 bits 1 hit · Pseudomonas · 177 bits WP_217662933.1 Arabidopsis lineage · rank 1 → 16 rank 1 · cellular organisms rank 2 · Eukaryota rank 3 · Viridiplantae ranks 4–11 · plants, angiosperms rank 12 · Pentapetalae ranks 13–15 · rosids, malvids rank 16 · Brassicaceae rule fires: age = rank 1 "cellular organisms" — from one hit fourteen intervening ranks, every one of them empty taxonomic representativeness 30% cutoff 22% output: "possible contamination or HGT" — no age assigned drop the entry → Brassicaceae, rank 16

The failure the gate does not catch

Representativeness protects against a spurious hit making a gene look too old. Nothing in the method protects against the opposite error — real orthologs that the search cannot see, making a gene look too young. That is the error that matters for de novo claims, and it has a systematic cause.

alignment length (aa) → bitscore 0 100 200 66 200 300 the tempting reading 74 hits at 116 bits · Brassicaceae — "real" 8 hits at 50–55 bits · everything else — "noise" 45% identity ortholog score ≈ 0.85 bits per residue 66 aa · ceiling ≈ 55 bits a perfect ortholog cannot score higher 8 hits, aligned over all 66 residues Nyssa sinensis 66/66 41% Vitis riparia 66/66 42% Bauhinia variegata 62/97 50% colinear · shared motifs · full length jackhmmer MSA: 185 taxa · 25 orders · only 18 Brassicaceae Fabales 39 · order 16 · Lamiales 15 · Rosales 15 · Solanales 13 · Asterales 12 · Ranunculales 5 · Proteales 1 the last two are basal eudicots — outside the rank genEra assigned short + fast-evolving + weakly expressed = both "de novo candidate" and "undetectable ortholog"
Direction of the bias

Homology detection failure does not add noise symmetrically. It moves genes down the ladder, never up. The concepts page covers the null model that puts a probability on it, and the synteny-based argument against it. Any statement of the form "N% of genes in this genome are lineage-specific" is therefore an upper bound, and the size of the overestimate scales with how short and fast-evolving the genes are.

Five ways a genEra run can be wrong while reporting success

All five were hit on this cluster across ~40 runs. In every case genEra printed Enjoy your results!!! and exited 0, and in four of the five it wrote one output row per query gene — so neither the exit code nor a line count catches any of them.

FailureWhat the output looks likeGuard
which leak100% of genes "Absent from the DIAMOND results"--cleanenv; assert type which is a binary
R_TempDir failure100% absent, on one bad node onlypass an explicit scratch TMPDIR into the container
FASTSTEP3R OOM100% absent; seff shows ~100% memory efficiencysize memory from the results-table size; gate on content
walltime truncationfewer rows than query sequences, summary still writtencompare row count against grep -c '^>'
wrong -t taxidoutput looks perfect; the ladder is the wrong speciesprint the resolved scientific name and read it

The last one is the most dangerous, because it produces a clean, plausible, fully populated result. In the campaign behind this site it happened twice: one species was run under the taxid of a different species in the same genus, and another was run as a hybrid when it is not one. Both were caught by an audit that compared the manifest against NCBI, not by anything in the pipeline.