Skip to content

Commit c429e09

Browse files
committed
Bio::EMBL: add #co / #contig (CO line)
CON-division entries (e.g. chromosome-level assemblies built by joining other entries, such as WGS scaffolds -- human GRCh38 chromosome 1, accession CM000663, is a real example) do not embed their sequence in an SQ/sequence record; instead they describe it as a join of other entries in a CO record. This line was previously unhandled, so such entries silently looked like empty (0 bp) sequences via sq/seq, with no way to retrieve the assembly instruction. Add #co (aliased as #contig, mirroring the recently added Bio::GenBank#contig), returning the CO content as a String. As with GenBank's #contig, sq/seq behavior for CO-only entries is left unchanged (they keep returning an empty sequence); #co/#contig is a complementary accessor, not a replacement.
1 parent f3c3215 commit c429e09

2 files changed

Lines changed: 70 additions & 0 deletions

File tree

lib/bio/db/embl/embl.rb

Lines changed: 16 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -462,6 +462,22 @@ def seq
462462
alias naseq seq
463463
alias ntseq seq
464464

465+
# CO -- Returns contents of the CO (contig/construction) line as a
466+
# String. Note that CON-division entries (e.g. chromosome-level
467+
# assemblies constructed by joining other entries, such as WGS
468+
# scaffolds) do not embed the sequence itself in an SQ/sequence
469+
# record, and instead describe it as a join of other entries in a
470+
# CO record; in this case, sq/seq return an empty sequence while co
471+
# returns the assembly instruction string (e.g.
472+
# "join(BX000000.1:1..1000,gap(100),BX000001.1:1..2000)"). Entries
473+
# without a CO line return an empty String.
474+
#
475+
# CO Line; contig/construction information (>=0 per entry)
476+
def co
477+
field_fetch('CO')
478+
end
479+
alias contig co
480+
465481
#--
466482
# // Line; termination line (end; 1/entry)
467483
#++

test/unit/bio/db/embl/test_embl_new_part.rb

Lines changed: 54 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -64,4 +64,58 @@ def test_pr_absent
6464
assert_equal([], @obj.pr)
6565
end
6666
end # class TestEMBL_PR_line_absent
67+
68+
class TestEMBL_CO_line < Test::Unit::TestCase
69+
def setup
70+
# CON-division entries (e.g. chromosome-level assemblies built
71+
# by joining WGS scaffolds, such as human GRCh38 chromosome 1,
72+
# accession CM000663) do not embed the sequence itself in an
73+
# SQ/sequence record; a CO line describes it as a join of other
74+
# entries instead.
75+
text = <<~THE_END_OF_THE_TEXT
76+
ID AB123456; SV 1; linear; genomic DNA; CON; PRO; 3100 BP.
77+
XX
78+
AC AB123456;
79+
XX
80+
CO join(BX000000.1:1..1000,gap(100),BX000001.1:1..2000)
81+
//
82+
THE_END_OF_THE_TEXT
83+
84+
@obj = Bio::EMBL.new(text)
85+
end
86+
87+
def test_co
88+
assert_equal('join(BX000000.1:1..1000,gap(100),BX000001.1:1..2000)',
89+
@obj.co)
90+
end
91+
92+
def test_contig
93+
assert_equal('join(BX000000.1:1..1000,gap(100),BX000001.1:1..2000)',
94+
@obj.contig)
95+
end
96+
97+
def test_seq_is_empty_when_only_co_line_present
98+
assert_equal('', @obj.seq.to_s)
99+
end
100+
end # class TestEMBL_CO_line
101+
102+
class TestEMBL_CO_line_absent < Test::Unit::TestCase
103+
def setup
104+
text = <<~THE_END_OF_THE_TEXT
105+
ID AB123456; SV 1; linear; genomic DNA; STD; PRO; 4 BP.
106+
XX
107+
AC AB123456;
108+
XX
109+
SQ Sequence 4 BP; 1 A; 1 C; 1 G; 1 T; 0 other;
110+
acgt 4
111+
//
112+
THE_END_OF_THE_TEXT
113+
114+
@obj = Bio::EMBL.new(text)
115+
end
116+
117+
def test_co_absent
118+
assert_equal('', @obj.co)
119+
end
120+
end # class TestEMBL_CO_line_absent
67121
end # module Bio

0 commit comments

Comments
 (0)