Revision history for Chem-Structure-Parser
0.038 2026-10-05 CDT
- structure_rmsd() keeps each row of its RMSD and count matrices as it
makes it, rather than fetching it back out of the matrix to fill it in.
Nothing it returns changes; it just no longer relies on a fetch it never
checked.
- The documentation now says what a residue's rsa is: its relative
solvent accessibility, its surface as a fraction of its maximum.
0.037 2026-10-01 CDT
- bzip2 files are read as gzip files are: a .bz2 is unpacked as it is
read, by every function that takes a file name -- structure_info() and
its views, structure_sequences() and structure_rmsd() -- in both formats,
with a file of several bzip2 streams, as pbzip2 writes, read to the end.
A .bz2 used to die saying it was not unpacked. A .Z still dies, now
saying that gzip and bzip2 are what is read instead.
- The minimum perl is now 5.10.1, up from 5.10.0: 5.10.1 is the first perl
whose core has IO::Uncompress::Bunzip2, which reads the .bz2 files.
0.036 2026-09-30 CDT
- A residue or serial number too big for its PDB columns is read as the
hybrid-36 number cctbx, phenix and gemmi write it as: residue A000 is
10000 and serial A0000 is 100000, as gemmi reads them. They used to be
no number at all, so every such residue in a chain was filed under the
same empty key -- three waters numbered A000, A001 and A002 came back as
one residue whose oxygen had three alternate conformers -- and a chain
past 9,999 residues read differently from its mmCIF.
- A residue with no number is no longer folded into a residue numbered 0
just before it, in either format. Both were read as residue 0, so the
second's atoms became conformers of the first's.
- A bad option to the feature calculations is refused before anything is
computed, rather than after the surface had already been written into
the structure.
- Reading is somewhat cheaper. The atom hashes are built with their keys
hashed once per file and at the size they end at, and an mmCIF header's
tables with their column names made once per table rather than once per
cell: under callgrind, 5% fewer instructions to read 4fqr and 9% fewer
to read the mmCIF of 2iw3. The feature walk reads each residue's
backbone and rings once where it read them three and four times, and
the interface reuses the surface's neighbour grid, for 4% fewer over
structure_features() on 1a4k. A file is read into a buffer of its own
size rather than one doubled up to it.
0.035 2026-09-27 CDT
- The interface has a pocket: the first partner's atoms within 6 A of the
second, which is the receptor's binding site around a peptide or ligand
in the order the partners are worked out. It comes with the pocket's
size, charge, Kyte-Doolittle hydropathy, residue-class fractions, surface
in the receptor alone and its hydrophobic share, mean relative surface
and buried fraction, and packing density. The definitions are the
pocket block of the pepPriML static feature script, and on 1cka all
sixteen of its figures agree. It is on by default; pocket => 0 turns it
off, and pocket_distance and packing_distance move its two cutoffs.
0.034 2026-09-26 CDT
- A free amino acid or nucleotide bound after a chain is told from a
modified residue capping it by whether it is bonded to the chain, which is
gemmi's test, and no longer by its residue number. The number let 1hsl's
free HIS 239 (one past a 238-residue chain) and 3lms's free GLY 501
(inside a chain numbered 4, 567, 1501, ...) into the sequence, and ten
PDBbind chains came back with a negative n_missing; none does now. It
also stops demoting the residues of a peptide written partly as HETATM:
cyclosporin's N-methyl amino acids are part of its sequence again.
- DBREF1/DBREF2 pairs are read. A cross-reference too long for DBREF's
columns -- a ten-character UniProt accession -- is written as the pair,
and 232 PDBbind entries that have nothing else had no dbref at all.
- A structure refined with SHELXL has an R-work and R-free. Its REMARK 3
spells them 'R VALUE (WORKING SET, NO CUTOFF)' and 'FREE R VALUE (NO
CUTOFF)', which were not recognised; the no-cutoff pair is the one the
archive's mmCIF files as the R values, and an mmCIF file whose
ls_R_factor_R_work is '?' now falls back to the same number in
_pdbx_refine.
- A REVDAT continuation line is part of the revision above it. Each one
used to come back as a revision of its own with no date, in 1,789
PDBbind entries; the records it lists are now added to its revision.
- mmCIF: a helix's id is the HELIX record's (AA1, or 1) rather than the
HELX_P1 the mmCIF file names its row; a sheet's strands carry their
n_strands and sense, from _struct_sheet and _struct_sheet_order, as a
SHEET record's do; a chain has the fragment of its entity; a synthetic
entity has the organism _pdbx_entity_src_syn names; and an entity from
several sources -- a chimera -- names each organism once, joined as
SOURCE joins them, rather than only the last row's.
- mmCIF: a chain id longer than seven characters is kept whole. It was cut
to seven, so two segments named SEG1PROA and SEG1PROB came back as one
chain.
- An option that needs a value and is given undef -- model => undef, probe
=> undef -- dies naming the function and the option. It died of 'Use of
uninitialized value' inside the module, which named neither. undef is
still off for a switch such as hydrogens, and still means not given for
chains and format.
- dssp => 1 with features => 0 computes the secondary structure on its own
and leaves it at {dssp}, and so does structure_info($file, 'dssp',
features => 0). The first returned nothing and said nothing; the second
died. dssp with atoms => 0, where there is nothing to compute it from,
dies.
- A coordinate field with anything after the number is not a number. A
coordinate too wide for its eight columns pushes its last digit into the
next field, and y came back 0 from '0-200.00'; it is now undef.
- The internal parse functions refuse a tied chains hash. One whose EXISTS
died took the parse down in the middle and leaked everything it had
built. structure_info() never passes one.
- The solvent-accessible surface is two to three times faster, and the
default read of a structure about twice as fast: 4fqr's from 2.47 s to
1.18 s, and 60 PDBbind entries from 10.77 s to 6.07 s. The sphere points
are tried in bins by direction, and a bin one neighbour covers whole is
passed over. Every surface is the same number to the last bit, which was
checked on 56 structures at three point counts and on the long-double and
quadmath perls.
- atoms => 0 reads in a third less memory and time. The parse still built
sixteen values for every atom that nobody would read, and now builds the
four the residues are worked out from: 2wy2's 64 models peak at 140 MB
rather than 408 MB, and 4fqr's read at 57 MB rather than 98 MB. What
comes back is the same structure.
- Reading with features => 0 is about a seventh faster: the atoms are filed
into their residues, alternate conformers and all, by the parse as it
builds them, rather than by a Perl loop over them afterwards. 60 PDBbind
entries went from 1.17 s to 1.00 s, and what comes back is the same.
- Edge-to-face pi stacking no longer depends on which ring comes first. The
intersect distance was measured from the first ring's centroid projected
onto the line where the two planes meet, which is mdtraj's arithmetic
asked one way round, and missed stacks whose second ring lies on that
line: 1b6c's PHE A99 is 0.66 A from it and its stack with TRP A59 was not
reported. It is now the nearer centroid's distance from the line, which
is what the documentation always said it was.
- Half-sphere exposure agrees with Biopython's HSExposureCB on every residue
of a hundred PDBbind entries but one, a neighbour 0.00005 degrees from the
dividing plane. Two rules were missing. A standard residue with no
peptide-linked standard neighbour -- the GLU of 1a08's bound ACE-FTY-GLU-
DIP -- is in none of CaPPBuilder's polypeptides, so it now gets no figure
and no longer adds one to the counts around it. And a side chain
deposited without its CB now gets no figure either, where it used to be
given glycine's virtual CB. Biopython takes the highest-occupancy
conformer, so altloc => 'highest' is the read that matches it.
- Deuterium is an element. A neutron structure's D atoms were counted in
n_no_element and given the 2.0 A default radius and no mass, which made an
N-D surface 30% larger than the N-H one; D now has hydrogen's radius and
mdtraj's deuterium mass, 2.0135532127.
- mmCIF: a quote its line never closes is read as part of a bare word, as
'CA, rather than as the start of a value that ran on through every line
below it until a stray quote turned up. The rows after one are no longer
lost into a single atom's name.
- mmCIF: a position modelled in two chemical states at once has one row in
_entity_poly_seq per state, and both were counted. 1ejg's SEQRES came
back 48 residues long from the mmCIF and 46 from the PDB file, with an
n_missing of 2 against 0; both now say 46 and 0.
- A gzip file of several members, which is what bgzip writes, is read to
the end. Only its first member was, so a structure split that way came
back short or with no atoms at all.
- A .bz2 or .Z file dies saying it is compressed with something this module
does not unpack. The format was worked out from the name with the suffix
taken off, and the compressed bytes then read as they stood, which came
back as a structure with no atoms and no error.
- structure_sequences($file) reads the file with atoms => 0 and
features => 0, as its sequences need neither. It had been computing the
surface and every other feature and throwing them away: over 233 files it
now takes 11.5 s where it took 64.
- On a perl whose NV is a double and whose compiler evaluates double
arithmetic in the x87's wider format -- gcc on 32-bit x86 -- the fast
number reader leaves fractions to strtod(). Its division was rounded
twice there, and 103 of 400,001 six-decimal strings came back one ulp away
from strtod()'s answer; no one- to five-decimal string tried was affected.
- The features include the interface between two partners, at
$info->{features}{interface}, and structure_interface() is new: the
residue pairs in contact across it and their classes, the surface the two
bury (apolar and polar), salt bridges, polar N/O/S pairs, cation-pi pairs,
bridging waters, the centre-of-mass and closest-atom distances, and the
B-factors of each side and of its interface. It is on by default, like
the rest, with the partners worked out -- the shortest polymer chain
against the others, or the largest ligand against a single chain -- or
named with partners => [ [...], [...] ] to structure_info(),
structure_features() or structure_interface(). interface => 0 turns it
off and sasa => 0 leaves out its surface. It makes the default read
6.4 times a bare one, where it was 4.7.
- The interface carries PRODIGY's non-interacting surface and binding
affinity estimate (prodigy-prot 2.4.0). Over 148 PDBbind entries the
contacts and their classes are PRODIGY's exactly, and the estimate is
within 0.22 kcal/mol of it; the difference is the two surfaces' radii.
- Each protein chain has a charge at pH 7 and an isoelectric point, which
are Biopython's to the bit, and each chain with secondary structure an
ss_fraction of H, E and C. The whole structure has the charges added up
and the fractions over all its chains.
0.033 2026-09-25 CDT
- structure_info($file, 'torsions') -- or 'torsion' -- returns the torsion
angles and nothing else: each chain's torsions hash, keyed by chain id, so
$t->{A}{phi} and $t->{B}{alpha} come from the same call when chain A is a
protein and chain B a nucleic acid. Each chain keeps only the angles its
residues have, and carries its residue_order, which the arrays run
parallel to and which is otherwise thrown away with the rest of the
structure.
- The chain id is column 22 of an ATOM or HETATM record and nothing else.
The reader used to fall back to column 21 when column 22 was blank, which
the format leaves empty and CHARMM and NAMD fill with the fourth letter of
a residue name: a TIP3 water with no chain id came back as chain '3', and
a POPC as chain 'C'. It never did the job it was there for either -- a
two-character chain id fills column 22 as well, so "AB" came back "B"
regardless. Column 22 alone is Biopython's rule. No entry in PDBbind
v2020 puts anything in column 21, so no deposited file reads differently.
- A formal charge written sign first in a PDB file, as 4byf and 4ui0 write
their O-1 oxygens, now reads as "1-", the format's own spelling and what
the mmCIF reader has always returned for the same atom. The two formats
gave the same entry two different charges.
- mmCIF: a pdbx_formal_charge too large to be a charge is no charge. An
integer as wide as an IV overflowed when its sign was taken off and came
back as "0-", and one too wide for an IV was cut to its first two
characters; both are now empty, as 10 always was.
- mmCIF: a text field in a file with DOS line ends reads the same as it does
from the same file with Unix ones. Every CR stayed in the value, one
before each line break and one at the end.
- A parse that dies while reading its options -- an option that is an object
whose numeric overload throws -- no longer leaks the hashes it had made
before it died, nor the file it had read.
- A structure with one coordinate that reads as nan no longer makes the
surface calculation quadratic. The neighbour grid took its box over every
atom, so one non-finite coordinate put the whole structure into one cell:
the same answer, the slow way -- 14.7 s against 2.0 s for 4fqr.
- Reading is faster and smaller. The seven short strings every atom carries
-- name, altloc, residue name, chain, insertion code, element, charge --
no longer have a buffer each, which was the largest cost in the parse:
structure_info() with features => 0 over 200 entries of PDBbind v2020 went
from 5.18 s to 4.38 s, and the peak of the default read of 4fqr, 90,792
atoms, from 234 MB to 217 MB. notes.txt has the measurements.
0.032 2026-09-21 CDT
- A path that is a directory is refused by the parse itself, on every system.
structure_info() has always said so, but the XS underneath it was left to
whatever the C library did with one, and that is three different things:
glibc opens it and fails the first read, Windows will not open it at all,
and NetBSD hands back the raw directory blocks, which read as a structure
with no atoms in it and nothing said so. Two 0.031 smokers -- Strawberry
5.42.0 and perl 5.42.3 on NetBSD 11 -- reported the two t/errors.t failures
this is the end of.
- The readers do not consult the locale. perl sets LC_CTYPE from the
caller's environment at startup, so <ctype.h> answered for it: in a Latin-1
locale a stray accented byte in the atom name column counted as a letter and
was guessed at as an element, and in a Turkish one an mmCIF file whose tags
are written in capitals would not have been recognised as mmCIF at all,
because tolower('I') there is not 'i'. Everything these readers classify --
an element symbol, a residue name, a CIF keyword -- is ASCII by definition,
and is now classified by perl's ASCII-only macros. The answer for a given
file no longer depends on the machine's locale, and t/parse.t reads a
structure again under every locale the machine has to say so.
- mmCIF: a type_symbol that does not spell an element is no longer believed.
The PDB reader has always refused columns 77-78 unless they are letters,
because an entry from before the element column existed keeps the entry id
there, and guessed the element from the atom name instead; a file converted
from one carries the same string in type_symbol, and a writer with no symbol
to write sometimes writes '' rather than '.'. Both now fall through to the
same guess, so the two formats agree about what the atom is made of.
- structure_info(file, features => { sasa => 0 }) is fatal rather than
ignored. The per-feature switches belong to structure_features(); spelled
into structure_info() the hash is simply a true value, so every feature was
computed -- the surface included, which is nearly all of the time a read
costs -- and the caller was told nothing.
- The solvent-accessible surface is about 9% faster, and every number it
returns is unchanged to the last bit. The kernel tests each of an atom's
thousand sphere points against the neighbours that might hide it and stops
at the first one that does, so it now tries the nearest neighbour first: a
nearer atom hides more of the sphere. Which neighbour is found first is all
that changes -- what is counted is the points no neighbour covers at all --
and over a 60-structure spread of PDBbind v2020 all 369,621 atom, residue
and chain surfaces compare byte for byte with 0.031's. The surface is
nearly all of what features cost, so a default structure_info() of a
directory comes out about 5% ahead.
- The whole-file read goes through PerlIO rather than stdio, which is the
layer perl's own opens use on every target, and a read that fails partway
through a plain file is reported rather than returned as a short one.
0.031 2026-09-19 CDT
- Fix for Windows: 0.03 did not compile there at all. A local variable in
contacts_find() was named near, and <windows.h>: which perl.h includes on
Win32 -- still defines near and far as empty macros, left over from 16-bit
segment addressing. The preprocessor deleted the name, the declaration
became a syntax error, and gcc rejected every later use of it. The
variable is renamed; nothing else about the contacts changes. http://matrix.cpantesters.org/?dist=Chem-Structure-Parser;os=mswin32;perl=5.42.2;reports=1
0.03 2026-09-18 CDT
- new: structure_rmsd(), how far apart two copies of a molecule are. Takes
file names, structures structure_info() already returned, or any mix of the
two, and compares every one of them with every other: two give the RMSD in
angstrom, more than two give the matrix. A file read with model => 'all'
counts as one structure per model, so an NMR ensemble is compared against
itself in one call --
my $r = structure_rmsd('2ll7.ent.pdb', model => 'all');
printf '%.2f', $r->{rmsd}[0][1]; # 3.61 A, models 1 and 2
Atoms are paired on the identity the file gives them -- chain, residue,
atom name -- so two models of an ensemble pair exactly and two files of the
same protein pair over whatever they share; nothing is aligned and nothing
is guessed. chains reads only some of them, chain_map says what a chain of
the later structures is called in the first, and match => 'order' pairs the
nth atom of each instead. select => 'heavy', 'backbone' or 'ca' narrows
what takes part, fit => 0 measures the two where they lie rather than
superposing them first, min_atoms is where "too few to mean anything" is
drawn, and detail => 1 returns the rotation and translation as well as the
number.
The superposition is Theobald's quaternion characteristic polynomial, after
qcprot.c and Biopython's Bio/PDB/qcprot.py. The deviation is not read off
the eigenvalue the way QCP's whole point is that it can be: that subtracts
two numbers which agree in as many figures as the two structures do, and two
structures being nearly the same is the ordinary case -- models 24 and 25 of
1JM4 have identical coordinates and gemmi 0.7.5, which takes that route,
answers 8.6e-07 A for them where this answers 0. Forming the rotation and
measuring with it costs one pass more over the paired atoms and has no
cancellation in it. Against gemmi's superpose_positions over 5,598
superpositions -- every pair of models of the 40 NMR entries in the first
two thousand files of PDBbind v2020 -- the largest relative difference is
5.65e-12, and against Biopython's SVDSuperimposer the two agree to every
figure either prints.
All of it is in C, for the reason the rest of the per-atom work is: over
2L6J, 20 models of 1,841 atoms and 190 pairs, 0.05 s against 1.19 s for the
same matrix written in Perl with the same arithmetic, and half the working
memory. benchmark.pl measures it and notes.txt records it.
- structure.info.pl gains --rmsd, the matrix of every structure named on the
command line against every other, and --select to say which atoms take
part. One NMR entry on its own compares its own models, because a file
with more than one counts as one structure per model.
- t/data gains ensemble.pdb and ensemble.cif, four models of chain B of 1K9R
as deposited: an ensemble whose models really differ, which nmr.pdb's three
do not -- they are one tripeptide moved a whole angstrom along x, and a
rigid move superposes onto itself exactly and says nothing about the fit.
t/rmsd.t is the new test file and carries gemmi's answer for every pair of
them, under each of the four selections, so that the comparison runs on a
machine with no Python on it.
- new: the rest of the physical properties, all computed by structure_info()
on the way past and reachable through $info->{features} or the functions
named below. Each is a translation of somebody else's implementation and is
tested against it rather than against itself.
- shape: the gyration tensor, its principal moments, and the asphericity,
acylindricity and relative shape anisotropy built from them, after
mdtraj's geometry/shape.py. Agrees with it to eight figures.
- interface: each chain's surface with the other chains taken away, and so
what the chains bury between them -- $f->{sasa}{buried}, and sasa_alone
and buried on every chain. One more surface calculation per chain, each
touching only that chain's atoms, so all of them together cost about what
the first one cost. interface => 0 turns it off.
- torsion angles: phi, psi, omega and chi1 to chi5 on every residue, in
degrees, after mdtraj's compute_phi and its relatives. Every one of them
is gathered onto its chain as well, as $chain->{torsions}{phi} and its
relatives: one array per angle, parallel to the chain's residue_order,
with an undef where a residue has no such angle. That is the form a
Ramachandran plot or a rotamer census wants, and it is the same option --
dihedrals => 0 leaves it off with the angles themselves.
- new: structure_contacts(), which residues touch which, after mdtraj's
compute_contacts() with its closest-heavy scheme; n_contacts on every
residue. Cross-chain and neighbouring pairs are reported as well as the
ones mdtraj's `all' considers.
- new: structure_hbonds(), the backbone hydrogen bonds by Kabsch and
Sander's electrostatic energy, after mdtraj's kabsch_sander(). The amide
hydrogen is placed rather than read, so it works on a crystal structure
that has none. Exact: the same bonds and the same energies to 4e-5
kcal/mol over 1A42, 1A22, 1AHW and 3AU6.
- secondary structure: ss and ss_simple on every residue, the Kabsch-Sander
dictionary. This is mdtraj's compute_dssp() letter for letter -- its
dssp.cpp (itself DSSP 2.2.0, ported by Robert T. McGibbon) transcribed
function for function, on a hydrogen bond table built to mdtraj's rules
rather than to this module's, computing in the float32 and the nanometre
mdtraj computes in because the assignment turns on two comparisons against
constants and one bond either way is worth several residues' letters.
t/features.t demands equality on every residue of every structure in
t/data rather than bounding a disagreement.
Measured over every tenth entry of PDBbind v2020 -- 1,011 of the 1,012
read, the other being a file mdtraj will not open -- the two give the same
letter for all 619,067 residues that have a backbone.
- new: structure_dssp(), the same assignment read the other way round --
chain id, then DSSP letter, then the positions in that chain's
residue_order of the residues that have it. That is the shape to take a
whole fold in at once, where ss on each residue is the shape to read while
walking them. structure_info($file, 'dssp') reads a file and hands back
just that hash, and structure_info($file, dssp => 1) leaves it at
$info->{dssp} alongside everything else.
A plain string in second place to structure_info() is new, and names a
view: the file is read, the view is taken out of it, and the rest is
thrown away. It cannot be confused with an option -- a file name followed
by an even number of arguments was never anything but a mistake before.
- half-sphere exposure: hse_up and hse_down on every standard amino acid
residue, after Biopython's HSExposureCB. Exact.
- nucleic acid torsions: alpha, beta, gamma, delta, epsilon, zeta and the
glycosidic chi on every nucleotide, plus nu, the five torsions of the
sugar ring, all in degrees. The atoms of each are the IUPAC-IUB (1983)
definitions; the angles agree with mdtraj's compute_dihedrals and gemmi's
calculate_dihedral over the same four atoms to 1.9e-4 degrees, which is
mdtraj's float32 coordinates, and to the rounding of the frozen file
against gemmi. alpha, epsilon and zeta span two residues and are not
measured across a chain break: the two must be joined by a phosphodiester
bond first, which is gemmi's are_connected() test at 2.4 A, and
phosphodiester_bond is the option that sets it.
- the sugar pucker: pucker_phase, pucker_amplitude and the envelope name
pucker on every nucleotide, by Altona and Sundaralingam's (1972)
pseudorotation equations, and glycosidic, which says whether the base is
anti or syn. This is the number that tells an A-form helix from a B-form
one, and it does: every ribose of the new rna.pdb is C3'-endo and every
deoxyribose of the new duplex.pdb is in the southern half of the cycle.
- the two numbers a nucleic acid sequence answers: gc_fraction and
purine_fraction on the structure and on every DNA or RNA chain, with
base_counts, n_gc and nucleotide_length beside them. Which letters count
is Bio.SeqUtils.gc_fraction()'s default, so an ambiguous base is in
neither half of the fraction rather than counted as something it is not.
- new: structure_base_pairs(), the Watson-Crick and wobble base pairs the
coordinates show, with the hydrogen bonds, the centroid separation, the
plane angle and the stagger each was found by, and base_pair on both of
the residues. Three pair types, numbered as Saenger (1984) numbers them:
19 for G-C, 20 for A-U and A-T, 28 for the G-U wobble and the G-T that is
the same pair in DNA.
This is the one property in the module with no reference implementation to
compare against -- mdtraj, gemmi and Biopython all stop short of base
pairing -- so it is measured against the annotation the wwPDB deposits
with the entry, which is 3DNA's. Over the forty archive entries that
carry an _ndb_struct_na_base_pair loop it finds all 1354 of the annotated
pairs that are between two unmodified bases and misses none. Eleven of
the eighteen it does not see have a modified base on one side and seven
are pairs with a crystallographic symmetry mate; two are found that the
annotation does not list, and both are argued in the documentation rather
than tuned away. base_pair_hbond and base_pair_stagger move the two
thresholds; base_pairs => 0 skips the work.
- new: structure_base_stacks(), how stacked every nearby pair of
nucleobases is, on the three variables the RNA literature scores stacking
with -- d0 between the two bases' centres of mass, the overlap angle
omega, and Xi between the two bases' normals -- and as one score from
-100% to 100% built out of them. base_stack goes onto both residues, as
base_pair does, and says which end of the pair each of them is, because
omega is measured from the 5' base and is not the same number read the
other way round.
The definition and every constant are Condon, D E; Kennedy, S D; Mort,
B C; Kierzek, R; Yildirim, I; Turner, D H (2015) "Stacking in RNA: NMR of
Four Tetramers Benchmark Molecular Dynamics", J Chem Theory Comput
11(6):2729-2742, section 2.4, and the implementation its numbers came out
of, which is the same authors' PDB_stacker. Figure 4 of the paper works
the three variables on residues 13 and 14 of 157D and reports 4.5 A,
40.7 degrees and 17.3 degrees; this answers 4.5286, 40.7402 and 17.3071,
and t/stacking.t holds it there against the new t/data/aform.pdb, which
is those residues as deposited. That one pair also settles the three
places where the paper's text and its authors' code disagree -- how
equation 10 is read, whether guanine's N2 is in its centre of mass, and
whether the distance knee is 3.5 A or 4 -- each of which is marked in
Parser.xs with the measurement that decides it.
Every base the module names is covered, DNA and RNA and the modified
spellings alike: a PSU is measured as a uridine and a 7MG as a guanosine,
through four sets of atoms serving the six letters A, C, G, I, T and U.
A base missing one of the atoms its entry names is in no pair rather than
measured wrong. base_stack_distance and base_stack_omega move the two
cutoffs; base_stacks => 0 skips the work.
- a torsion angle and a Kabsch-Sander amide hydrogen are no longer measured
across a chain break. mdtraj takes the residue before this one to be
whichever came before it in the file and computes a phi, or places a
hydrogen, across whatever gap is there; here the two must be peptide-bonded
first, at the 1.8 A Biopython's PPBuilder uses. Four collinear atoms get no
torsion either, rather than the zero atan2(0, 0) returns.
- t/data gains ss.pdb, fold.pdb and sheet.pdb, with their mmCIF twins: five
cysteines of 1AHW, two of them bonded and one not; sixty consecutive residues
of 1A42 that hold an alpha helix, a 3-10 helix and the turns and bends around
them; and sixty of 1H1V that hold, when read on their own, every one of the
eight DSSP letters -- the two strands of a ladder, an isolated bridge, and
the pi helix, which is the one letter a structure has to be looked for to
find. Without the second the hydrogen bonds would have had almost nothing in
t/data to be tested on, and without the third the beta half of the secondary
structure -- the bridges, the ladders and the bulges between them, which is
most of the code -- would have had nothing at all.
- and rna.pdb, duplex.pdb and wobble.pdb, likewise with their twins: six
nucleotides of 1MSY's GUAA tetraloop, which is where the riboses are, four
base pairs of 1BNA, the Drew-Dickerson dodecamer, which is where B-form DNA
is, and twelve nucleotides of 1MSY's lower helix, which is where the wobble
is. bases.pdb is a single DNA strand and has no O2' in it, so without the
first two the nucleic acid torsions would have been tested on half of what
they say; the third holds all three of the pair types base_pairs() looks for
and one that it has to leave alone.
- structure.info.pl gains --ss, one row per disulfide with what the file
declares beside it.
- structure_info() computes the physical properties on the way past and
leaves them in $info->{features}, so most callers never name
structure_features() at all. This is not free: reading goes from about
214,000 atoms a second to about 19,500, eleven times the cost, nearly all of
it the solvent-accessible surface at 960 sphere points per atom. With
sasa => 0 everything else together is 2.7x the read, and interface => 0
alone takes the whole thing from 10.9x to 7.1x.
features => 0 is the way back to the old speed,
and is what to reach for when reading a directory for its headers or its
sequences. atoms => 0 turns them off on its own rather than dying, so
structure_sequences($file, atoms => 0) still works.
- with no options, structure_features(), structure_sasa(),
structure_pi_stacking() and structure_disulfides() hand back what the read
already computed rather than walking the structure a second time. Name any
option and they compute again with it in force.
- new: structure_disulfides(), and $features->{disulfides} -- the disulfide
bonds the coordinates show, as against $info->{ssbond}, which is what the
depositor said. Each bond is also written onto both of its residues, as
{disulfide}. The rule is mdtraj's: a cysteine with an SG and no HG, paired
with another under 3.0 A, so a cysteine whose thiol hydrogen was modelled is
reduced and holds no bond. disulfide_distance moves the cutoff.
The two answers agree on every entry with disulfides in a 60-entry spread of
PDBbind, and their bond lengths agree to the two decimals SSBOND is written
in. Over a wider sweep they agree on 21 of 22: 1A4K is a Fab that is in the
file twice, and its SSBOND records cover one copy while its coordinates show
both.
The comparison is against mdtraj's rule and not against mdtraj's function.
Topology.create_disulfide_bonds() is handed the angstrom positions its PDB
reader parsed and compares them with 0.3, meaning nanometres, so it tests
for an SG-SG separation under 0.3 A and finds nothing on any file -- on
1A4K, 8 pairs under 0.3 nm and 0 under 0.3 A. The SG-SG bonds that do turn
up in its topology come from the file's own CONECT records.
- fixed: a structure whose atom records have unreadable coordinate columns was
mistaken for one read with atoms => 0 and died where it should have carried
on. The two are told apart by the atom names the residues list, not by the
positions that came back.
- new: structure_features(), the physical properties of a structure in one
call -- solvent-accessible surface, radius of gyration, molecular weight,
centroid and centre of mass, mean Kyte-Doolittle hydropathy, the aromatic
fraction of the sequence, and the stacked pairs of aromatic rings. One XS
walk of the structure serves all of them, because each has to touch every
atom and flattening a hash of chains of residues of atoms into coordinate
arrays is the expensive part.
- new: structure_sasa(), the solvent-accessible surface on its own, and
structure_pi_stacking(), the stacked ring pairs on their own.
- the per-atom, per-residue and per-chain results are written into the
structure beside the atom, residue and chain they belong to:
$info->{chains}{A}{sasa}, {residues}{54}{sasa}, {residues}{54}{rsa} --
the residue's surface as a fraction of the most it could have -- and
{residues}{54}{atoms}{CZ}{sasa}. store => 0 leaves the structure alone
and returns the totals only.
- the arithmetic is a translation of somebody else's, and is tested against
it rather than against itself: the Shrake-Rupley surface, the van der
Waals radii, the atomic masses, the ring geometry and the radius of
gyration are mdtraj 1.11's, and the hydropathy and aromaticity scales are
Biopython's. t/features.t compares every atom, every residue and every
stacked pair against what mdtraj answered for the structures in t/data,
frozen in t/data/features.txt, and re-runs mdtraj where it is installed so
the frozen answer cannot go stale. Run in float64, mdtraj's own kernel
and this one agree on every atom to nine digits.
- the neighbour search is a uniform grid rather than the all-pairs scan
mdtraj uses, which is what makes the surface of a large entry finish: the
cost is flat at about 24 microseconds an atom from 6,000 atoms to 16,000,
where an all-pairs scan grows with the square. Both find the same
neighbours, which is what the comparison above establishes.
- two departures from mdtraj, both deliberate and both documented where they
are made. Its face-to-face centroid cutoff of 5.5 is nanometres in a
function whose other distances are nanometres -- 55 A, far enough that any
two rings in a small protein would qualify -- where ProLIF, which the
geometry comes from, has 5.5 A; face_distance defaults to 5.5 A and
face_distance => 55 gets mdtraj's number. And compute_rg(traj, masses)
weights distances by mass but measures them from the geometric centroid;
rg_mass measures from the centre of mass, while rg keeps mdtraj's
equal-weight default, where the two centres coincide.
- a ligand contributes no aromatic rings. Finding one means perceiving its
bonds, which this module does not do -- it reads no CONECT record and
guesses no bond -- so structure_pi_stacking() says so rather than quietly
finding nothing.
- structure.info.pl gains --features, one tab separated row of properties per
structure, and --stacks, one row per stacked pair.
- removed: pdb_info() and cif_info(). Each was structure_info() with
format => 'pdb' or format => 'mmcif' filled in, and nothing else -- no
option, no error and no answer of their own -- so the two names bought a
reader nothing but the question of whether they differed. Detection costs
nothing to skip: it reads the file only when the name says nothing at all.
Call structure_info($file, format => 'pdb') where you called pdb_info($file),
and structure_info($file, format => 'mmcif') where you called cif_info($file).
- the release notes are no longer generated from README.md. They lived in
three places at once -- the README, the module's POD and this file -- and
the POD copy was a hundred and fifty lines of history in the middle of the
reference documentation. This file is now written by hand and md2pod.pl
only checks it.
- one fewer prerequisite: autodie is no longer required. It was there to
check five calls -- two opens, a read and two closes -- and those are now
checked where they are made, dying in autodie's own words, so a failure
says exactly what it used to say. Nothing a caller can see changes.
- fixed: a truncated or corrupt .gz was read as a short file and not as the
error it is. IO::Uncompress::Gunzip's read() returns a negative number for
one, where the loop reading the archive took anything not positive for the
end of the stream: half of a gzipped entry came back as a structure with no
atoms in it and nothing said so. It now dies, naming the file and what the
decompressor said was wrong with it. A directory of downloads is where such
a file is found, and a silent half-read of one is the worst answer available.
- fixed: an integer field too big for an IV wrapped instead of being declined.
No PDB record can reach it -- a serial number has five columns and a residue
number four -- but an mmCIF value is free-form, and twenty digits of
_atom_site.id came back as 200376420520689663 and twenty-one digits of
auth_seq_id as a negative residue number. A number that does not fit is now
undef, which is what '*****' in an overflowed PDB column already gave, and
IV_MAX and IV_MIN still read back exactly.
- fixed: a numeric field of digits and dots that is not a number took the whole
read down. The pattern behind CRYST1's cell, a HELIX length and every mmCIF
number matched '1.2.3' and then added it to zero, which under
warnings FATAL => 'all' is fatal rather than merely wrong. Fields of that
shape are in the archive -- 5m04 writes its pH as '5.4.-5.8' -- and one now
reads as the undef that says the file gave no number, leaving the fields
beside it alone.
- fixed: model => 2 of a file with no MODEL records came back with no atoms and
no chains and said nothing. Asking for a model that is not there falls back
to the first one, but the fall-back read the list of model numbers the file
declared, and a crystal structure declares none: its one model is model 1.
It now falls back for those too, which is most of the archive. A structure
emptied by something else -- a chains list that matches nothing, a file with
no coordinates -- is untouched and is not read a second time.
- structure_contacts() no longer takes time proportional to the residues times
the atoms. Each residue scanned the whole heavy-atom array to find its own,
though the atoms of a residue are one contiguous run of it; walking that run
instead took the contacts of 4fqr -- 11,406 residues, 88,416 heavy atoms --
from 2.25 seconds to 0.15, and the whole feature walk of that entry from 5.98
to 3.87, the rest of which is the solvent-accessible surface. It finds the
same pairs: t/features.t compares them with mdtraj's, as it did before, and
they are the same pairs an all-against-all walk of the fixtures finds.
- the parse builds a third fewer SVs. A residue's identity -- its name, its
chain, its number, its insertion code, whether it is a HETATM and which model
it is in -- was built once per atom and read once per residue. It is now
built where it is read, which for the shape structure_info() asks for is six
strings and integers per residue instead of six per atom: over 200 entries of
PDBbind v2020 the coordinate section went from 1.95 to 1.59 seconds, and
reading the largest of them (4fqr, 90,792 atoms) from 168 MB of peak memory
to 141 and from 0.28 seconds to 0.21. A caller of the low-level parse that
asks for columns still gets one entry per atom, because there the columns are
the answer; t/parse.t reads both shapes of four fixtures and asserts they say
the same thing.
- the surface of each chain alone -- the interface half of the properties, and
about a third of what they cost -- no longer computes every atom twice. An
atom whose every neighbour is in its own chain has the same surface with the
other chains there and taken away, bit for bit, so only the atoms that do
have a neighbour across a chain boundary are computed again. Over 20 entries
of PDBbind v2020 the interface pass went from 1.85 seconds to 0.38, and the
whole default feature walk of 4fqr from 4.38 seconds to 2.61. It is the same
answer: t/features.t still compares each chain's isolated surface with
mdtraj's, and t/properties.t compares it with reading that chain on its own,
which agreed to the last bit on all 131 chains of the 32 multi-chain entries
of a PDBbind v2020 sample as well as on the fixtures.
- the sphere-point loop inside the surface reads its neighbours out of four
arrays of its own rather than through an index into the structure's, with
each neighbour's radius squared once when it is found rather than once per
point. Same comparisons, same answers, about 9% off the surface.
- fixed: missing_residues could name a residue the chain has. A polymer
numbered out of order has jumps that pass over each other, and the list was
the jumps' contents added up rather than the set of numbers the chain does
not have -- so a residue modelled out of line could appear in its own chain's
missing list, and a number two jumps passed over could appear twice. It is
now a set, ascending. No chain of the 1,000 entries of PDBbind v2020 this
was measured over is affected; a chain numbered out of order is.
- structure_features() no longer trusts a chain's n_atoms to size its walk. It
is an upper bound for every structure this module reads, and a structure a
caller assembled by hand can say less than it holds, or not say at all --
where the walk quietly stopped at the count instead of at the atoms.
- the test suite exercises a good deal more of the module: every name in the
residue table rather than a sample of them (150 of them, checked against
gemmi 0.7.5 where gemmi's own table carries the name), the mmCIF categories
and syntax the fixtures have no use for, the header records an entry has one
of, the low-level parse in both of its shapes, and the error paths of the XS
entry points. Statement coverage of the Perl half goes from 96.5% to 99.6%
and of the XS from 94.3% to 96.6%, and branch coverage from 82.2% to 88.9%
and from 83.1% to 85.6%.
- t/leaks.t no longer fails under Devel::Cover. Test::LeakTrace does not work
with it and has a guard that says so, but Devel::Cover 1.52 installs no
runops routine for the guard to notice, so every statement it had
instrumented was counted as a leak and `cover -t' failed 65 of that file's 85
tests on Devel::Cover's own counters. The file now makes the skip itself,
running each block anyway so that a coverage run still covers it, and the
reference-cycle checks at the foot of it still run for real.
- fixed: installing the module no longer installs Test::More and Test::Exception
with it. dist.ini listed both in a [Prereqs] block with no phase, which is
the runtime phase, so the released META.json asked every installer for the
two modules only t/*.t loads -- both releases so far, and not the checked-in
Makefile.PL, which has had them in TEST_REQUIRES all along.
They are in [Prereqs / TestRequires] now. Carp is gone from Makefile.PL's
PREREQ_PM for the same reason in reverse: it was listed and loaded nowhere.
- SECURITY.md says where to take a report that goes unanswered -- the CPAN
Security Group, <cpan-security@security.metacpan.org>, who can triage it and
who are the better first stop for anybody who would rather not write to an
individual maintainer -- and which perls the support statement covers, which
is every perl the module installs on and not merely the newest.
0.02 2026-08-22 CDT
- Fix for Windows: https://www.cpantesters.org/cpan/report/960c8de0-8001-1014-8897-e6a949c55ebf
0.01 2026-08-21 CDT
- initial version: reads PDB into a hash of hashes, single-letter
sequences, residue types. Release notes from here on are the
maintainer's to write.
- reads mmCIF/PDBx (.cif, .mmcif, .pdbx) as well as PDB. structure_info()
works out which, and returns the same hash of hashes either way: the same
chains, residues, atoms, sequences and counts. Tested by reading fixture
pairs both ways and comparing them with is_deeply, and by converting real
PDB entries to mmCIF and asserting that nothing changes.
- the mmCIF reader uses the auth_* identifiers, so a chain read from a .cif
has the same name and residue numbering as the same chain read from a
.pdb, and converts the values the two formats spell differently (a formal
charge of -1 reads back as '1-', and a charge of 0 stays '0', which is
not the same answer as a blank charge field).
- new: cif_info(), the mmCIF counterpart of pdb_info().
- new: $info->{chains}{$chain}{missing_residues}, the residue numbers the
chain's gaps step over, in ascending order. Every chain has the list;
a chain with no gaps has it empty. It counts numbers rather than
residues, so a chain numbered by homology to a reference protein --
chymotrypsin numbering and the conventions like it, about one chain in
twenty-five -- reports the numbers its scheme skips on purpose along
with the ones that went unmodelled. n_missing is the count to trust
when a file has SEQRES and the two disagree.
- gaps no longer reads a change of numbering scheme as a gap. An antibody
numbered by the Kabat scheme runs 27, 1027, 2027, 28, where the
thousands are insertions after 27 and not a 999-residue hole, and 1a4k
was reading as a 214-residue light chain missing five thousand
residues. A chain can only be missing as many residues as the span
from its first polymer residue to its last leaves room for, counting
insertion codes as the one number they share, and a jump wider than
that is no longer counted. Nothing is read from SEQRES to decide it,
so a chain answers the same whether it came from a PDB file or an
mmCIF one and whether or not the headers were parsed.
- formats() now reports mmcif as supported.
- new: $info->{stats}{total_atoms}, every ATOM and HETATM record the file
has, every model and before the model selection or the hydrogens,
waters, hetatm and chains options threw anything away. n_atoms is what
came back and this is what there was, so total_atoms == n_atoms +
n_skipped whatever the options were set to, and a structure read out of
a 64-model ensemble can say that its 6,432 atoms are one model of
411,648 rather than the whole file. Both readers count it the same way.
- new: is_single_ion($chain), or is_single_ion($info, $chain), true when a
chain holds exactly one residue. An ion given a chain of its own is a
chain with no sequence to read, and a structure with a dozen of them has
more of those chains than polymer ones, so the loop that puts them aside
is worth not writing by hand. In XS, and single counts residues in the
chain: not atoms in the residue, so a sulphate and a perchlorate answer
the same, and not the residue's type, which comes off a table of names
that cannot be complete -- SO4 is on the module's ION list and BF4 is
not, and that is a fact about the list. The residue is not asked what it
is, so a chain of one sugar or one water reads true as well; the residue
says which it is, in its own type. A chain of two zincs is not one, and
neither is a protein chain with a zinc numbered into it. Handing it the
whole structure with no chain id, or a residue, is fatal: all three are
hash references and a false answer would be taken at face value.
- fixed: SEQRES was read to the end of the line rather than to column 70.
An entry deposited before about 1996 keeps its id and a line number in
columns 73-80 of every record, and those became two more residues per
SEQRES line: pdb1gdr's 140-residue chain read as 162 residues with an X
every thirteenth place, which is a wrong sequence rather than a missing
one. No file in the remediated archive is affected; the ones the archive
still distributes as deposited are.
- fixed: the element columns are no longer believed when they do not spell
an element. The same files put part of the entry id in columns 77-78, so
every atom of pdb1gdr read as element '1' -- which also stopped
hydrogens => 0 from finding hydrogens, since the element is what says
which atoms those are. A field that is not letters falls back to the
atom name, which the module already knows how to read. Likewise the
charge columns: a charge is a digit and a sign, and 'DR' is not one, so
it reads as the empty string a blank field would have given.
- fixed: a HELIX length that is not a number now reads as empty rather
than as the text that was in columns 72-76.
- fixed: the text records -- TITLE, COMPND, SOURCE, KEYWDS, AUTHOR, EXPDTA,
JRNL -- are cut at column 72 when columns 73-80 hold nothing but the
entry id and a line number. Text that is not the entry id is left alone,
so a title that really does run to column 80 is not truncated.
- fixed: a free-text COMPND or SOURCE is no longer thrown away. A file
older than the MOL_ID convention writes 'COMPND GAMMA DELTA RESOLVASE'
and names no chains, so the entry is the one molecule and every chain in
it gets it; $info->{compound}{1}{free_text} says the record was read that
way rather than parsed into tokens.
- fixed: an atom's altlocs list was missing the conformer that supplied the
coordinates when that record had no altloc letter. disordered.pdb writes
ARG 27's CZ once with a blank altloc and once as B; the list held only
the B, so the occupancies of an atom summed to 0.5 and a caller writing
the conformers back out wrote one of two.
- resolution now falls back to REMARK 3's RESOLUTION RANGE HIGH when there
is no REMARK 2. A file written by a refinement program rather than by
the archive often has the whole of REMARK 3 and no REMARK 2 at all, and
it is the same number the mmCIF reader already takes from
_refine.ls_d_res_high, so the two formats answer alike. REMARK 2 still
wins where there is one, and a BIN RESOLUTION RANGE HIGH is never it.
- new: t/foreign.t, the cases that gemmi's and Biopython's own test
directories know about -- an atom whose first record has no altloc
letter, a coordinate line that stops early, the same record written
twice, a MODEL with no serial number, a serial number that spills out of
its columns, one residue modelled in two chemical states in mmCIF, the
element rules with no element columns, a resolution that is only in
REMARK 3, CRLF line endings, and CIFs that are not structures.
- new: t/data/pdb1gdr.ent, a 1993 entry that keeps its id in columns 73-80.
- new: t/oracle.t, which compares every atom of model 1 against gemmi --
chain, number, insertion code, name, altloc and coordinates, as a
multiset -- over t/data and a spread of STRUCTURE_INFO_TEST_DIR and
STRUCTURE_INFO_TEST_CIF_DIR. t/real.t checks the C against a second
reader in Perl, which cannot catch a column both of them read wrongly
because one person wrote both. It skips unless python3 can import gemmi.
Over 655 real structures the two agree atom for atom except where a
residue is modelled in two chemical states at once, which is one residue
here and two there.
- new: $info->{chains}{$chain}{elements}, how many atoms of each element the
chain holds. Same shape as $info->{stats}{elements}, which is the whole
structure; both count coordinate records, as the n_atoms beside them
does, so both add up to it. Gathered in the parse, once per residue for
the lookup and once per atom for the increment, so it costs the same as
the whole-structure tally already did. With model => 'all' each model's
chains carry their own.
- element symbols now read back as IUPAC writes them: Zn, not ZN. Columns
77-78 of a PDB record are capitals, an mmCIF type_symbol is capitals as
often as not, and guess_element() uppercases what it takes from the atom
name, so a zinc arrived as ZN by all three roads. The correction runs
once, on the symbol, where it is settled, so an atom's element, the chain
tally and the structure tally cannot disagree. Only the 118 named
elements are corrected; a field that spells no element is left as the
file wrote it, so XX stays XX rather than becoming a plausible Xx.
Incompatible: $atom->{element}, the element column of the low-level
parse, and the keys of $info->{stats}{elements} all change spelling for
the two-letter elements.
- the element counts are unsigned integers rather than whatever sv_inc()
left behind. They are counted up from nothing and never down.
Keyboard Shortcuts
Global
s
Focus search bar
?
Bring up this help dialog
GitHub
gp
Go to pull requests
gi
Go to GitHub issues (only if GitHub is preferred repository)