diff --git a/RcppTskit/inst/examples/create_test.trees.R b/RcppTskit/inst/examples/create_test.trees.R
index 02cbfa2..213da2d 100644
--- a/RcppTskit/inst/examples/create_test.trees.R
+++ b/RcppTskit/inst/examples/create_test.trees.R
@@ -20,6 +20,8 @@ builtins <- import_builtins()
 msprime <- import("msprime")
 tskit <- import("tskit")
 
+# -----------------------------------------------------------------------------
+
 # Generate a tree sequence for testing
 ts <- msprime$sim_ancestry(
   samples = 80,
@@ -85,6 +87,8 @@ length(ts$tables$individuals$metadata) # 0
 ts$dump("inst/examples/test.trees")
 # ts <- tskit$load("inst/examples/test.trees")
 
+# -----------------------------------------------------------------------------
+
 # Create a second tree sequence with metadata in some tables
 # basic_schema <- tskit$MetadataSchema("{'codec': 'json'}")
 # Can't get this to work via reticulate :(
@@ -104,3 +108,34 @@ ts$metadata_schema # {"codec":"json"}
 ts$tables$individuals$metadata # R vector
 builtins$type(ts$tables$individuals$metadata) # numpy.ndarray
 length(ts$tables$individuals$metadata) # 21
+
+# -----------------------------------------------------------------------------
+
+# Another example with a reference sequence
+
+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
+ts$reference_sequence # NULL
+
+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()
+iterate(ali, function(i) cat(i, "\n"))
+# iterate(ali, function(i) print(i)) # produces no output
+ali_vec <- iterate(ts$alignments())
+print(ali_vec)
+
+ts$dump("RcppTskit/inst/examples/test_with_ref_seq.trees")
+
+# -----------------------------------------------------------------------------
diff --git a/RcppTskit/inst/examples/create_test.trees.py b/RcppTskit/inst/examples/create_test.trees.py
index 5f94cdf..4c00823 100644
--- a/RcppTskit/inst/examples/create_test.trees.py
+++ b/RcppTskit/inst/examples/create_test.trees.py
@@ -2,6 +2,8 @@
 import tskit
 import os
 
+# -----------------------------------------------------------------------------
+
 # Generate a tree sequence for testing
 ts = msprime.sim_ancestry(
     samples=80, sequence_length=1e4, recombination_rate=1e-4, random_seed=42
@@ -55,11 +57,13 @@
 ts.tables.individuals.metadata.shape  # (0,)
 
 os.getcwd()
-os.chdir("RcppTskit")
-ts.dump("inst/examples/test.trees")
-# ts = tskit.load("test.trees")
+ts.dump("RcppTskit/inst/examples/test.trees")
+# ts = tskit.load("RcppTskit/inst/examples/test.trees")
+
+# -----------------------------------------------------------------------------
 
 # Create a second tree sequence with metadata in some tables
+# ts = tskit.load("RcppTskit/inst/examples/test.trees")
 ts2_tables = ts.dump_tables()
 len(ts2_tables.metadata)
 ts2_tables.metadata = tskit.pack_bytes('{"seed": 42, "note": "ts2"}')
@@ -107,7 +111,7 @@
 len(ts.tables.individuals.metadata)  # 21
 ts.tables.individuals.metadata.shape  # (21,)
 
-ts.dump("inst/examples/test2.trees")
+ts.dump("RcppTskit/inst/examples/test2.trees")
 
 tables = ts.dump_tables()
 tables.metadata_schema = tskit.MetadataSchema(None)
@@ -116,3 +120,26 @@
 ts = tables.tree_sequence()
 ts.metadata
 len(ts.metadata)  # 23
+
+# -----------------------------------------------------------------------------
+
+# Another example with a reference sequence
+
+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()
+ts.reference_sequence
+
+tables = ts.dump_tables()
+tables.reference_sequence.data = "ATCGAATTCG"
+ts = tables.tree_sequence()
+ts.has_reference_sequence()
+ts.reference_sequence
+
+ali = ts.alignments()
+for i in ali:
+    print(i)
+
+ts.dump("RcppTskit/inst/examples/test_with_ref_seq.trees")
+
+# -----------------------------------------------------------------------------
diff --git a/RcppTskit/inst/examples/test_with_ref_seq.trees b/RcppTskit/inst/examples/test_with_ref_seq.trees
new file mode 100644
index 0000000..3a94227
Binary files /dev/null and b/RcppTskit/inst/examples/test_with_ref_seq.trees differ
diff --git a/RcppTskit/tests/testthat/test_load_summary_and_dump.R b/RcppTskit/tests/testthat/test_load_summary_and_dump.R
index 40bac75..6881186 100644
--- a/RcppTskit/tests/testthat/test_load_summary_and_dump.R
+++ b/RcppTskit/tests/testthat/test_load_summary_and_dump.R
@@ -7,33 +7,57 @@ test_that("ts/tc_load(), ts/tc_summary*(), and ts/tc_dump(x) work", {
   expect_error(ts_load("nonexistent_ts"))
   ts_file <- system.file("examples/test.trees", package = "RcppTskit")
 
+  expect_error(
+    ts_ptr_load(ts_file, options = bitwShiftL(1L, 30)),
+    regexp = "ts_ptr_load only supports load options"
+    # TSK_LOAD_SKIP_TABLES (1 << 0) and TSK_LOAD_SKIP_REFERENCE_SEQUENCE (1 << 1)
+  )
+
+  # TSK_LOAD_SKIP_TABLES (1 << 0)
   expect_error(
     ts_load(ts_file, skip_tables = "y"),
     regexp = "skip_tables must be TRUE/FALSE!"
   )
-  expect_no_error(tc_load(ts_file, skip_tables = TRUE))
+  expect_no_error(ts_ptr_load(ts_file, options = bitwShiftL(1L, 0L)))
+  expect_no_error(ts_load(ts_file, skip_tables = TRUE))
+  check_empty_tables_ptr <- function(ts) {
+    # jarl-ignore implicit_assignment: it's just a test
+    tmp <- capture.output(p <- ts_ptr_print(ts))
+    expect_true(all(p$tables$number == 0))
+  }
   check_empty_tables <- function(ts) {
     # jarl-ignore implicit_assignment: it's just a test
     tmp <- capture.output(p <- ts$print())
     expect_true(all(p$tables$number == 0))
   }
+  ts_ptr <- ts_ptr_load(ts_file, options = bitwShiftL(1L, 0L))
+  check_empty_tables_ptr(ts_ptr)
   ts <- ts_load(ts_file, skip_tables = TRUE)
   check_empty_tables(ts)
-  ts <- TableCollection$new(file = ts_file, skip_tables = TRUE)
+  ts <- TreeSequence$new(file = ts_file, skip_tables = TRUE)
   check_empty_tables(ts)
 
+  # TSK_LOAD_SKIP_REFERENCE_SEQUENCE (1 << 1)
   expect_error(
     ts_load(ts_file, skip_reference_sequence = 1L),
     regexp = "skip_reference_sequence must be TRUE/FALSE!"
   )
+  expect_no_error(ts_ptr_load(ts_file, options = bitwShiftL(1L, 1L)))
   expect_no_error(ts_load(ts_file, skip_reference_sequence = TRUE))
 
-  expect_error(
-    ts_ptr_load(ts_file, options = bitwShiftL(1L, 30)),
-    regexp = "ts_ptr_load only supports load options"
-    # TSK_LOAD_SKIP_TABLES (1 << 0) and TSK_LOAD_SKIP_REFERENCE_SEQUENCE (1 << 1)
+  ts_with_ref_seq_file <- system.file(
+    "examples/test_with_ref_seq.trees",
+    package = "RcppTskit"
   )
-
+  expect_no_error(ts_ptr_load(ts_with_ref_seq_file))
+  expect_no_error(ts_ptr_load(
+    ts_with_ref_seq_file,
+    options = bitwShiftL(1L, 1L)
+  ))
+  expect_no_error(ts_load(ts_with_ref_seq_file))
+  expect_no_error(ts_load(ts_with_ref_seq_file, skip_reference_sequence = TRUE))
+
+  # For tests below
   ts_ptr <- ts_ptr_load(ts_file)
   ts <- ts_load(ts_file)
 
@@ -45,33 +69,57 @@ test_that("ts/tc_load(), ts/tc_summary*(), and ts/tc_dump(x) work", {
   expect_error(tc_load("nonexistent_ts"))
   ts_file <- system.file("examples/test.trees", package = "RcppTskit")
 
+  expect_error(
+    tc_ptr_load(ts_file, options = bitwShiftL(1L, 30)),
+    regexp = "tc_ptr_load only supports load options"
+    # TSK_LOAD_SKIP_TABLES (1 << 0) and TSK_LOAD_SKIP_REFERENCE_SEQUENCE (1 << 1)
+  )
+
+  # TSK_LOAD_SKIP_TABLES (1 << 0)
   expect_error(
     tc_load(ts_file, skip_tables = "y"),
     regexp = "skip_tables must be TRUE/FALSE!"
   )
+  expect_no_error(tc_ptr_load(ts_file, options = bitwShiftL(1L, 0L)))
   expect_no_error(tc_load(ts_file, skip_tables = TRUE))
+  check_empty_tables_ptr <- function(tc) {
+    # jarl-ignore implicit_assignment: it's just a test
+    tmp <- capture.output(p <- tc_ptr_print(tc))
+    expect_true(all(p$tables$number == 0))
+  }
   check_empty_tables <- function(tc) {
-    # jarl-ignore implicit_assignment:  it's just a test
+    # jarl-ignore implicit_assignment: it's just a test
     tmp <- capture.output(p <- tc$print())
     expect_true(all(p$tables$number == 0))
   }
+  tc_ptr <- tc_ptr_load(ts_file, options = bitwShiftL(1L, 0L))
+  check_empty_tables_ptr(tc_ptr)
   tc <- tc_load(ts_file, skip_tables = TRUE)
   check_empty_tables(tc)
   tc <- TableCollection$new(file = ts_file, skip_tables = TRUE)
   check_empty_tables(tc)
 
+  # TSK_LOAD_SKIP_REFERENCE_SEQUENCE (1 << 1)
   expect_error(
     tc_load(ts_file, skip_reference_sequence = 1L),
     regexp = "skip_reference_sequence must be TRUE/FALSE!"
   )
+  expect_no_error(tc_ptr_load(ts_file, options = bitwShiftL(1L, 1L)))
   expect_no_error(tc_load(ts_file, skip_reference_sequence = TRUE))
 
-  expect_error(
-    tc_ptr_load(ts_file, options = bitwShiftL(1L, 30)),
-    regexp = "tc_ptr_load only supports load options"
-    # TSK_LOAD_SKIP_TABLES (1 << 0) and TSK_LOAD_SKIP_REFERENCE_SEQUENCE (1 << 1)
+  ts_with_ref_seq_file <- system.file(
+    "examples/test_with_ref_seq.trees",
+    package = "RcppTskit"
   )
-
+  expect_no_error(tc_ptr_load(ts_with_ref_seq_file))
+  expect_no_error(tc_ptr_load(
+    ts_with_ref_seq_file,
+    options = bitwShiftL(1L, 1L)
+  ))
+  expect_no_error(tc_load(ts_with_ref_seq_file))
+  expect_no_error(tc_load(ts_with_ref_seq_file, skip_reference_sequence = TRUE))
+
+  # For tests below
   tc_ptr <- tc_ptr_load(ts_file)
   tc <- tc_load(ts_file)