Test if skip_ref_seq option works properly

[gregorgorjanc]

Will need a tree sequence with a ref sequence first! How do we get one?

@LynxJinyangii could you help me out here?

See the examples I have created at https://github.com/HighlanderLab/RcppTskit/tree/main/RcppTskit/inst/examples.

Perhaps best if we take the example.trees and simplify it down to something even smaller (take the exmaples code from

RcppTskit/RcppTskit/R/RcppTskit.R

Line 579 in 9af8e17

#' tmp <- tc_py$simplify(samples = c(0L, 1L, 2L, 3L))

) and we add ref seq to it. I am not clear how this is done! See https://tskit.dev/tskit/docs/stable/data-model.html#reference-sequence, https://tskit.dev/tskit/docs/stable/python-api.html#tskit.TreeSequence.alignments - this one has this code to get random sequence reference_sequence=tskit.random_nucleotides(ts.sequence_length), but we need to store it into ts or tc object - how do we do this?)

• Add ref seq into ts and create example_with_ref_seq.trees
• Develop tests

[LynxJinyangii]

ts = msprime.sim_ancestry(
    samples=3, ploidy=2, sequence_length=10, random_seed=2)
ts = msprime.sim_mutations(ts, rate=0.1, random_seed=2)
ts.has_reference_sequence()

False

tables = ts.dump_tables()
tables.reference_sequence.data = "ATCGAATTCG"
ts = tables.tree_sequence()
ts.has_reference_sequence()

True

ts.reference_sequence

ReferenceSequence({'metadata_schema': '', 'metadata': b'', 'data': 'ATCGAATTCG', 'url': ''})

ali = ts.alignments()
for i in ali:
    print(i)

ATACCATTCA
ATACGACTCA
ATACCATTCC
ATATCATTCA
ATGCGACTCA
ATGCGACTCA

Something could be interesting FYI (although I am not sure if we really need this feature for now...):
tskit-dev/tskit#146 (comment)
tskit-dev/tskit@43921cc

[gregorgorjanc]

Thanks @LynxJinyangii this is nice simple example!

[gregorgorjanc]

Have added the above example and unit test to ensure C++ and R arguments work. However, atm there is no way to test if we read in the data properly since we don't have accessor function to pull ref seq into R from the C/C++ structs. This is getting a bit annoying, which is why I am looking at #49.

Will leave this issue open, but rename it to be clearer.

• So, the next action is to check if the ref seq is correctly read in, but we probably want some work done or design decisions taken on Providing zero-copy read access to the table columns and tree arrays? #49 first (to avoid writting copy-getters for each property/element of ts&tc structs