diff --git a/RcppTskit/NEWS.md b/RcppTskit/NEWS.md
index 5ff9c90..652c623 100644
--- a/RcppTskit/NEWS.md
+++ b/RcppTskit/NEWS.md
@@ -50,6 +50,7 @@ and releases adhere to [Semantic Versioning](https://semver.org/spec/v2.0.0.html
 
 - Turn vignette URL as hyperlinks and similar cosmetics.
 - State mirroring of the `R/Python` APIs and `C++/C` APIs across the package.
+- Update tskit C to 1.3.1
 - TODO
 
 ## [0.2.0] - 2026-02-22
diff --git a/RcppTskit/inst/include/tskit/VERSION.txt b/RcppTskit/inst/include/tskit/VERSION.txt
index 589268e..6261a05 100644
--- a/RcppTskit/inst/include/tskit/VERSION.txt
+++ b/RcppTskit/inst/include/tskit/VERSION.txt
@@ -1 +1 @@
-1.3.0
\ No newline at end of file
+1.3.1
\ No newline at end of file
diff --git a/RcppTskit/inst/include/tskit/tskit/VERSION_kastore.txt b/RcppTskit/inst/include/tskit/tskit/VERSION_kastore.txt
index 3e3c2f1..eca07e4 100644
--- a/RcppTskit/inst/include/tskit/tskit/VERSION_kastore.txt
+++ b/RcppTskit/inst/include/tskit/tskit/VERSION_kastore.txt
@@ -1 +1 @@
-2.1.1
+2.1.2
diff --git a/RcppTskit/inst/include/tskit/tskit/core.h b/RcppTskit/inst/include/tskit/tskit/core.h
index 115ef6c..9512a90 100644
--- a/RcppTskit/inst/include/tskit/tskit/core.h
+++ b/RcppTskit/inst/include/tskit/tskit/core.h
@@ -42,7 +42,7 @@ extern "C" {
 
 #ifdef __GNUC__
 #define TSK_WARN_UNUSED __attribute__((warn_unused_result))
-#define TSK_UNUSED(x) TSK_UNUSED_##x __attribute__((__unused__))
+#define TSK_UNUSED(x)   TSK_UNUSED_##x __attribute__((__unused__))
 #else
 #define TSK_WARN_UNUSED
 #define TSK_UNUSED(x) TSK_UNUSED_##x
@@ -93,11 +93,11 @@ missing data.
  * on the thread above.
  */
 typedef int64_t tsk_id_t;
-#define TSK_MAX_ID INT64_MAX - 1
+#define TSK_MAX_ID          INT64_MAX - 1
 #define TSK_ID_STORAGE_TYPE KAS_INT64
 #else
 typedef int32_t tsk_id_t;
-#define TSK_MAX_ID INT32_MAX - 1
+#define TSK_MAX_ID          INT32_MAX - 1
 #define TSK_ID_STORAGE_TYPE KAS_INT32
 #endif
 
@@ -109,7 +109,7 @@ The ``tsk_size_t`` type is an unsigned integer used for any size or count value.
 @endrst
 */
 typedef uint64_t tsk_size_t;
-#define TSK_MAX_SIZE UINT64_MAX
+#define TSK_MAX_SIZE          UINT64_MAX
 #define TSK_SIZE_STORAGE_TYPE KAS_UINT64
 
 /**
@@ -152,7 +152,7 @@ to the API or ABI are introduced, i.e., the addition of a new function.
 The library patch version. Incremented when any changes not relevant to the
 to the API or ABI are introduced, i.e., internal refactors of bugfixes.
 */
-#define TSK_VERSION_PATCH   0
+#define TSK_VERSION_PATCH   1
 /** @} */
 
 /*
diff --git a/RcppTskit/inst/include/tskit/tskit/kastore.h b/RcppTskit/inst/include/tskit/tskit/kastore.h
index fbd48ce..3acc01a 100644
--- a/RcppTskit/inst/include/tskit/tskit/kastore.h
+++ b/RcppTskit/inst/include/tskit/tskit/kastore.h
@@ -157,7 +157,7 @@ to the API or ABI are introduced, i.e., the addition of a new function.
 The library patch version. Incremented when any changes not relevant to the
 to the API or ABI are introduced, i.e., internal refactors of bugfixes.
 */
-#define KAS_VERSION_PATCH   1
+#define KAS_VERSION_PATCH   2
 /** @} */
 
 #define KAS_HEADER_SIZE             64
diff --git a/RcppTskit/inst/include/tskit/tskit/tables.h b/RcppTskit/inst/include/tskit/tskit/tables.h
index 9523ee1..9618e4c 100644
--- a/RcppTskit/inst/include/tskit/tskit/tables.h
+++ b/RcppTskit/inst/include/tskit/tskit/tables.h
@@ -4607,7 +4607,7 @@ tsk_id_t tsk_table_collection_check_integrity(
 /* Undocumented methods */
 
 /* Flags for ibd_segments */
-#define TSK_IBD_STORE_PAIRS (1 << 0)
+#define TSK_IBD_STORE_PAIRS    (1 << 0)
 #define TSK_IBD_STORE_SEGMENTS (1 << 1)
 
 /* TODO be systematic about where "result" should be in the params
diff --git a/RcppTskit/inst/include/tskit/tskit/trees.h b/RcppTskit/inst/include/tskit/tskit/trees.h
index 2e495a4..84480ed 100644
--- a/RcppTskit/inst/include/tskit/tskit/trees.h
+++ b/RcppTskit/inst/include/tskit/tskit/trees.h
@@ -663,7 +663,7 @@ sequence.
 @endrst
 
 @param self A pointer to a tsk_treeseq_t object.
-@return Returns a pointer to the time units.
+@return Returns a pointer to the null-terminated file uuid.
 */
 const char *tsk_treeseq_get_file_uuid(const tsk_treeseq_t *self);
 
diff --git a/RcppTskit/src/tskit/VERSION.txt b/RcppTskit/src/tskit/VERSION.txt
index 589268e..6261a05 100644
--- a/RcppTskit/src/tskit/VERSION.txt
+++ b/RcppTskit/src/tskit/VERSION.txt
@@ -1 +1 @@
-1.3.0
\ No newline at end of file
+1.3.1
\ No newline at end of file
diff --git a/RcppTskit/src/tskit/VERSION_kastore.txt b/RcppTskit/src/tskit/VERSION_kastore.txt
index 3e3c2f1..eca07e4 100644
--- a/RcppTskit/src/tskit/VERSION_kastore.txt
+++ b/RcppTskit/src/tskit/VERSION_kastore.txt
@@ -1 +1 @@
-2.1.1
+2.1.2
diff --git a/RcppTskit/src/tskit/core.c b/RcppTskit/src/tskit/core.c
index 07a85b2..1f4f17d 100644
--- a/RcppTskit/src/tskit/core.c
+++ b/RcppTskit/src/tskit/core.c
@@ -883,8 +883,8 @@ tsk_search_sorted(const double *restrict array, tsk_size_t size, double value)
             upper = mid;
         }
     }
-    offset = (int64_t)(array[lower] < value);
-    return (tsk_size_t)(lower + offset);
+    offset = (int64_t) (array[lower] < value);
+    return (tsk_size_t) (lower + offset);
 }
 
 /* Rounds the specified double to the closest multiple of 10**-num_digits. If
diff --git a/RcppTskit/src/tskit/haplotype_matching.c b/RcppTskit/src/tskit/haplotype_matching.c
index a8acd20..d68776a 100644
--- a/RcppTskit/src/tskit/haplotype_matching.c
+++ b/RcppTskit/src/tskit/haplotype_matching.c
@@ -572,7 +572,7 @@ static inline bool
 element_in(
     const uint64_t *restrict A, tsk_id_t u, const tsk_id_t state, tsk_size_t num_words)
 {
-    tsk_size_t index = ((tsk_size_t) u) * num_words + (tsk_size_t)(state / 64);
+    tsk_size_t index = ((tsk_size_t) u) * num_words + (tsk_size_t) (state / 64);
     return (A[index] & (1ULL << (state % 64))) != 0;
 }
 
@@ -580,7 +580,7 @@ static inline void
 set_optimal_value(
     uint64_t *restrict A, tsk_id_t u, const tsk_size_t num_words, tsk_id_t state)
 {
-    tsk_size_t index = ((tsk_size_t) u) * num_words + (tsk_size_t)(state / 64);
+    tsk_size_t index = ((tsk_size_t) u) * num_words + (tsk_size_t) (state / 64);
     tsk_bug_assert(((tsk_size_t) state) / 64 < num_words);
     A[index] |= 1ULL << (state % 64);
 }
@@ -815,7 +815,7 @@ tsk_ls_hmm_redistribute_transitions(tsk_ls_hmm_t *self)
      * iterating over the roots. See the existing parsimony implementations
      * for an example. */
     for (root = tsk_tree_get_left_root(&self->tree); root != TSK_NULL;
-         root = right_sib[root]) {
+        root = right_sib[root]) {
         stack[0].tree_node = root;
         stack[0].old_state = T_old[T_index[root]].value_index;
         stack[0].new_state
@@ -966,7 +966,7 @@ tsk_ls_hmm_run_forward(tsk_ls_hmm_t *self, int32_t *haplotype)
     }
 
     for (t_ret = tsk_tree_first(&self->tree); t_ret == TSK_TREE_OK;
-         t_ret = tsk_tree_next(&self->tree)) {
+        t_ret = tsk_tree_next(&self->tree)) {
         ret = tsk_ls_hmm_update_tree(self, TSK_DIR_FORWARD);
         if (ret != 0) {
             goto out;
@@ -1166,7 +1166,7 @@ tsk_ls_hmm_run_backward(
     }
 
     for (t_ret = tsk_tree_last(&self->tree); t_ret == TSK_TREE_OK;
-         t_ret = tsk_tree_prev(&self->tree)) {
+        t_ret = tsk_tree_prev(&self->tree)) {
         ret = tsk_ls_hmm_update_tree(self, TSK_DIR_REVERSE);
         if (ret != 0) {
             goto out;
@@ -1492,7 +1492,7 @@ tsk_compressed_matrix_decode(tsk_compressed_matrix_t *self, double *values)
     }
 
     for (t_ret = tsk_tree_first(&tree); t_ret == TSK_TREE_OK;
-         t_ret = tsk_tree_next(&tree)) {
+        t_ret = tsk_tree_next(&tree)) {
         ret = tsk_tree_get_sites(&tree, &sites, &num_tree_sites);
         if (ret != 0) {
             goto out;
diff --git a/RcppTskit/src/tskit/tables.c b/RcppTskit/src/tskit/tables.c
index 8b92be4..d673619 100644
--- a/RcppTskit/src/tskit/tables.c
+++ b/RcppTskit/src/tskit/tables.c
@@ -100,7 +100,7 @@ check_offset_overflow(tsk_size_t current_size, tsk_size_t additional_elements)
            || current_size > (max_val - additional_elements);
 }
 
-#define TSK_NUM_ROWS_UNSET ((tsk_size_t) -1)
+#define TSK_NUM_ROWS_UNSET   ((tsk_size_t) - 1)
 #define TSK_MAX_COL_NAME_LEN 64
 
 static int
@@ -8255,8 +8255,8 @@ pair_to_integer(tsk_id_t a, tsk_id_t b, tsk_size_t N)
 static inline void
 integer_to_pair(int64_t index, tsk_size_t N, tsk_id_t *a, tsk_id_t *b)
 {
-    *a = (tsk_id_t)(index / (int64_t) N);
-    *b = (tsk_id_t)(index % (int64_t) N);
+    *a = (tsk_id_t) (index / (int64_t) N);
+    *b = (tsk_id_t) (index % (int64_t) N);
 }
 
 static int64_t
@@ -8997,7 +8997,7 @@ simplifier_check_state(simplifier_t *self)
     for (j = 0; j < self->input_tables.nodes.num_rows; j++) {
         last_position = -1;
         for (list_node = self->node_mutation_list_map_head[j]; list_node != NULL;
-             list_node = list_node->next) {
+            list_node = list_node->next) {
             tsk_bug_assert(
                 self->input_tables.mutations.node[list_node->mutation] == (tsk_id_t) j);
             site = self->input_tables.mutations.site[list_node->mutation];
@@ -9028,7 +9028,7 @@ simplifier_check_state(simplifier_t *self)
         child = self->buffered_children[j];
         tsk_bug_assert(self->child_edge_map_head[child] != NULL);
         for (int_list = self->child_edge_map_head[child]; int_list != NULL;
-             int_list = int_list->next) {
+            int_list = int_list->next) {
             tsk_bug_assert(int_list->left < int_list->right);
             if (int_list->next != NULL) {
                 tsk_bug_assert(int_list->right < int_list->next->left);
@@ -9108,7 +9108,7 @@ simplifier_print_state(simplifier_t *self, FILE *out)
         child = self->buffered_children[j];
         fprintf(out, "%lld -> ", (long long) j);
         for (int_list = self->child_edge_map_head[child]; int_list != NULL;
-             int_list = int_list->next) {
+            int_list = int_list->next) {
             fprintf(out, "(%f, %f), ", int_list->left, int_list->right);
         }
         fprintf(out, "\n");
@@ -9123,7 +9123,7 @@ simplifier_print_state(simplifier_t *self, FILE *out)
         if (self->node_mutation_list_map_head[j] != NULL) {
             fprintf(out, "%lld\t-> [", (long long) j);
             for (list_node = self->node_mutation_list_map_head[j]; list_node != NULL;
-                 list_node = list_node->next) {
+                list_node = list_node->next) {
                 fprintf(out, "%lld,", (long long) list_node->mutation);
             }
             fprintf(out, "]\n");
@@ -10187,7 +10187,7 @@ simplifier_insert_input_roots(simplifier_t *self)
     const double *node_time = self->tables->nodes.time;
 
     for (input_id = 0; input_id < (tsk_id_t) self->input_tables.nodes.num_rows;
-         input_id++) {
+        input_id++) {
         x = self->ancestor_map_head[input_id];
         if (x != NULL) {
             output_id = self->node_id_map[input_id];
@@ -10736,11 +10736,11 @@ tsk_table_collection_check_individual_integrity(
 
     for (j = 0; j < (tsk_size_t) num_individuals; j++) {
         for (k = individuals.parents_offset[j]; k < individuals.parents_offset[j + 1];
-             k++) {
+            k++) {
             /* Check parent references are valid */
             if (individuals.parents[k] != TSK_NULL
                 && (individuals.parents[k] < 0
-                       || individuals.parents[k] >= num_individuals)) {
+                    || individuals.parents[k] >= num_individuals)) {
                 ret = tsk_trace_error(TSK_ERR_INDIVIDUAL_OUT_OF_BOUNDS);
                 goto out;
             }
@@ -11250,20 +11250,20 @@ tsk_table_collection_equals(const tsk_table_collection_t *self,
     if (!(options & TSK_CMP_IGNORE_TABLES)) {
         ret = ret
               && tsk_individual_table_equals(
-                     &self->individuals, &other->individuals, options)
+                  &self->individuals, &other->individuals, options)
               && tsk_node_table_equals(&self->nodes, &other->nodes, options)
               && tsk_edge_table_equals(&self->edges, &other->edges, options)
               && tsk_migration_table_equals(
-                     &self->migrations, &other->migrations, options)
+                  &self->migrations, &other->migrations, options)
               && tsk_site_table_equals(&self->sites, &other->sites, options)
               && tsk_mutation_table_equals(&self->mutations, &other->mutations, options)
               && tsk_population_table_equals(
-                     &self->populations, &other->populations, options);
+                  &self->populations, &other->populations, options);
         /* TSK_CMP_IGNORE_TABLES implies TSK_CMP_IGNORE_PROVENANCE */
         if (!(options & TSK_CMP_IGNORE_PROVENANCE)) {
             ret = ret
                   && tsk_provenance_table_equals(
-                         &self->provenances, &other->provenances, options);
+                      &self->provenances, &other->provenances, options);
         }
     }
     /* TSK_CMP_IGNORE_TS_METADATA is implied by TSK_CMP_IGNORE_METADATA */
@@ -11273,19 +11273,19 @@ tsk_table_collection_equals(const tsk_table_collection_t *self,
     if (!(options & TSK_CMP_IGNORE_TS_METADATA)) {
         ret = ret
               && (self->metadata_length == other->metadata_length
-                     && self->metadata_schema_length == other->metadata_schema_length
-                     && tsk_memcmp(self->metadata, other->metadata,
-                            self->metadata_length * sizeof(char))
-                            == 0
-                     && tsk_memcmp(self->metadata_schema, other->metadata_schema,
-                            self->metadata_schema_length * sizeof(char))
-                            == 0);
+                  && self->metadata_schema_length == other->metadata_schema_length
+                  && tsk_memcmp(self->metadata, other->metadata,
+                         self->metadata_length * sizeof(char))
+                         == 0
+                  && tsk_memcmp(self->metadata_schema, other->metadata_schema,
+                         self->metadata_schema_length * sizeof(char))
+                         == 0);
     }
 
     if (!(options & TSK_CMP_IGNORE_REFERENCE_SEQUENCE)) {
         ret = ret
               && tsk_reference_sequence_equals(
-                     &self->reference_sequence, &other->reference_sequence, options);
+                  &self->reference_sequence, &other->reference_sequence, options);
     }
     return ret;
 }
@@ -12057,7 +12057,7 @@ tsk_table_collection_dumpf(
         { "format/name", (const void *) &TSK_FILE_FORMAT_NAME,
             TSK_FILE_FORMAT_NAME_LENGTH, KAS_INT8 },
         { "format/version",
-            (const void *) &(uint32_t[]){
+            (const void *) &(uint32_t[]) {
                 TSK_FILE_FORMAT_VERSION_MAJOR, TSK_FILE_FORMAT_VERSION_MINOR },
             2, KAS_UINT32 },
         { "sequence_length", (const void *) &self->sequence_length, 1, KAS_FLOAT64 },
@@ -13302,7 +13302,7 @@ tsk_table_collection_union(tsk_table_collection_t *self,
     /* Now we know the full individual map we can remap the parents of the new
      * individuals*/
     for (k = (tsk_id_t) self->individuals.parents_offset[num_individuals_self];
-         k < (tsk_id_t) self->individuals.parents_length; k++) {
+        k < (tsk_id_t) self->individuals.parents_length; k++) {
         if (self->individuals.parents[k] != TSK_NULL) {
             self->individuals.parents[k] = individual_map[self->individuals.parents[k]];
         }
diff --git a/RcppTskit/src/tskit/trees.c b/RcppTskit/src/tskit/trees.c
index fdc043e..1aa06e5 100644
--- a/RcppTskit/src/tskit/trees.c
+++ b/RcppTskit/src/tskit/trees.c
@@ -228,7 +228,7 @@ tsk_treeseq_init_individuals(tsk_treeseq_t *self)
         if (ind != TSK_NULL) {
             node_array = self->individual_nodes[ind];
             tsk_bug_assert(node_array - self->individual_nodes_mem
-                           < (tsk_id_t)(total_node_refs - node_count[ind]));
+                           < (tsk_id_t) (total_node_refs - node_count[ind]));
             node_array[node_count[ind]] = node;
             node_count[ind] += 1;
         }
@@ -537,7 +537,7 @@ tsk_treeseq_init(
 
     if (tsk_treeseq_get_time_units_length(self) == strlen(TSK_TIME_UNITS_UNCALIBRATED)
         && !strncmp(tsk_treeseq_get_time_units(self), TSK_TIME_UNITS_UNCALIBRATED,
-               strlen(TSK_TIME_UNITS_UNCALIBRATED))) {
+            strlen(TSK_TIME_UNITS_UNCALIBRATED))) {
         self->time_uncalibrated = true;
     }
 out:
@@ -856,7 +856,7 @@ tsk_treeseq_get_individuals_time(const tsk_treeseq_t *self, double *output)
 
 /* Stats functions */
 
-#define GET_2D_ROW(array, row_len, row) (array + (((size_t)(row_len)) * (size_t)(row)))
+#define GET_2D_ROW(array, row_len, row) (array + (((size_t) (row_len)) * (size_t) (row)))
 
 static inline double *
 GET_3D_ROW(double *base, tsk_size_t num_nodes, tsk_size_t output_dim,
@@ -899,7 +899,7 @@ tsk_treeseq_genealogical_nearest_neighbours(const tsk_treeseq_t *self,
     /* TODO It's probably not worth bothering with the int16_t here. */
     int16_t k, focal_reference_set;
     /* We use the K'th element of the array for the total. */
-    const int16_t K = (int16_t)(num_reference_sets + 1);
+    const int16_t K = (int16_t) (num_reference_sets + 1);
     tsk_size_t num_nodes = self->tables->nodes.num_rows;
     const tsk_id_t num_edges = (tsk_id_t) self->tables->edges.num_rows;
     const tsk_id_t *restrict I = self->tables->indexes.edge_insertion_order;
@@ -1082,7 +1082,7 @@ tsk_treeseq_mean_descendants(const tsk_treeseq_t *self,
     tsk_size_t j;
     int32_t k;
     /* We use the K'th element of the array for the total. */
-    const int32_t K = (int32_t)(num_reference_sets + 1);
+    const int32_t K = (int32_t) (num_reference_sets + 1);
     tsk_size_t num_nodes = self->tables->nodes.num_rows;
     const tsk_id_t num_edges = (tsk_id_t) self->tables->edges.num_rows;
     const tsk_id_t *restrict I = self->tables->indexes.edge_insertion_order;
@@ -1742,7 +1742,7 @@ tsk_treeseq_site_general_stat(const tsk_treeseq_t *self, tsk_size_t state_dim,
 
         /* Update the sites */
         for (tree_site = 0; tree_site < self->tree_sites_length[tree_index];
-             tree_site++) {
+            tree_site++) {
             site = self->tree_sites[tree_index] + tree_site;
             ret = compute_general_stat_site_result(site, state, state_dim, result_dim, f,
                 f_params, total_weight, polarised, site_result);
@@ -2882,7 +2882,7 @@ positions_to_tree_indexes(const tsk_treeseq_t *ts, const double *positions,
         }
         (*tree_indexes)[i] = tree_index;
     }
-    tsk_bug_assert(tree_index <= (tsk_id_t)(num_trees - 1));
+    tsk_bug_assert(tree_index <= (tsk_id_t) (num_trees - 1));
 
 out:
     return ret;
@@ -2896,7 +2896,7 @@ get_index_counts(
     tsk_id_t index = indexes[0];
     tsk_size_t count, i;
     tsk_size_t *counts = tsk_calloc(
-        (tsk_size_t)(indexes[num_indexes ? num_indexes - 1 : 0] - indexes[0] + 1),
+        (tsk_size_t) (indexes[num_indexes ? num_indexes - 1 : 0] - indexes[0] + 1),
         sizeof(*counts));
     if (counts == NULL) {
         ret = tsk_trace_error(TSK_ERR_NO_MEMORY);
@@ -3315,7 +3315,7 @@ tsk_treeseq_two_branch_count_stat(const tsk_treeseq_t *self, tsk_size_t state_di
         }
         col = 0;
         for (c = 0; c < (col_indexes[n_cols ? n_cols - 1 : 0] - col_indexes[0] + 1);
-             c++) {
+            c++) {
             ret = advance_collect_edges(&r_state, (tsk_id_t) c + col_indexes[0]);
             if (ret != 0) {
                 goto out;
@@ -3621,7 +3621,7 @@ tsk_treeseq_site_allele_frequency_spectrum(const tsk_treeseq_t *self,
 
         /* Update the sites */
         for (tree_site = 0; tree_site < self->tree_sites_length[tree_index];
-             tree_site++) {
+            tree_site++) {
             site = self->tree_sites[tree_index] + tree_site;
             while (windows[window_index + 1] <= site->position) {
                 window_index++;
@@ -4566,8 +4566,8 @@ D2_unbiased_summary_func(tsk_size_t state_dim, const double *state,
         double w_ab = n - (w_AB + w_Ab + w_aB);
         result[j] = (1 / (n * (n - 1) * (n - 2) * (n - 3)))
                     * ((w_aB * w_aB * (w_Ab - 1) * w_Ab)
-                          + ((w_ab - 1) * w_ab * (w_AB - 1) * w_AB)
-                          - (w_aB * w_Ab * (w_Ab + (2 * w_ab * w_AB) - 1)));
+                        + ((w_ab - 1) * w_ab * (w_AB - 1) * w_AB)
+                        - (w_aB * w_Ab * (w_Ab + (2 * w_ab * w_AB) - 1)));
     }
     return 0;
 }
@@ -4604,8 +4604,8 @@ Dz_unbiased_summary_func(tsk_size_t state_dim, const double *state,
         result[j] = (1 / (n * (n - 1) * (n - 2) * (n - 3)))
                     * ((((w_AB * w_ab) - (w_Ab * w_aB)) * (w_aB + w_ab - w_AB - w_Ab)
                            * (w_Ab + w_ab - w_AB - w_aB))
-                          - ((w_AB * w_ab) * (w_AB + w_ab - w_Ab - w_aB - 2))
-                          - ((w_Ab * w_aB) * (w_Ab + w_aB - w_AB - w_ab - 2)));
+                        - ((w_AB * w_ab) * (w_AB + w_ab - w_Ab - w_aB - 2))
+                        - ((w_Ab * w_aB) * (w_Ab + w_aB - w_AB - w_ab - 2)));
     }
     return 0;
 }
@@ -4639,11 +4639,10 @@ pi2_unbiased_summary_func(tsk_size_t state_dim, const double *state,
         double w_Ab = state_row[1];
         double w_aB = state_row[2];
         double w_ab = n - (w_AB + w_Ab + w_aB);
-        result[j]
-            = (1 / (n * (n - 1) * (n - 2) * (n - 3)))
-              * (((w_AB + w_Ab) * (w_aB + w_ab) * (w_AB + w_aB) * (w_Ab + w_ab))
-                    - ((w_AB * w_ab) * (w_AB + w_ab + (3 * w_Ab) + (3 * w_aB) - 1))
-                    - ((w_Ab * w_aB) * (w_Ab + w_aB + (3 * w_AB) + (3 * w_ab) - 1)));
+        result[j] = (1 / (n * (n - 1) * (n - 2) * (n - 3)))
+                    * (((w_AB + w_Ab) * (w_aB + w_ab) * (w_AB + w_aB) * (w_Ab + w_ab))
+                        - ((w_AB * w_ab) * (w_AB + w_ab + (3 * w_Ab) + (3 * w_aB) - 1))
+                        - ((w_Ab * w_aB) * (w_Ab + w_aB + (3 * w_AB) + (3 * w_ab) - 1)));
     }
     return 0;
 }
@@ -7784,7 +7783,7 @@ tsk_tree_map_mutations(tsk_tree_t *self, int32_t *genotypes,
             ret = tsk_trace_error(TSK_ERR_BAD_ANCESTRAL_STATE);
             goto out;
         } else if (ancestral_state >= num_alleles) {
-            num_alleles = (int32_t)(ancestral_state + 1);
+            num_alleles = (int32_t) (ancestral_state + 1);
         }
     }
 
@@ -7882,8 +7881,8 @@ kc_vectors_alloc(kc_vectors *self, tsk_id_t n)
 
     self->n = n;
     self->N = (n * (n - 1)) / 2;
-    self->m = tsk_calloc((size_t)(self->N + self->n), sizeof(*self->m));
-    self->M = tsk_calloc((size_t)(self->N + self->n), sizeof(*self->M));
+    self->m = tsk_calloc((size_t) (self->N + self->n), sizeof(*self->m));
+    self->M = tsk_calloc((size_t) (self->N + self->n), sizeof(*self->M));
     if (self->m == NULL || self->M == NULL) {
         ret = tsk_trace_error(TSK_ERR_NO_MEMORY);
         goto out;
@@ -10739,14 +10738,14 @@ tsk_matvec_calculator_run(tsk_matvec_calculator_t *self)
     while (m < self->num_windows) {
         if (valid && self->position == tree_pos.interval.left) {
             for (k = (tsk_size_t) tree_pos.out.start;
-                 k != (tsk_size_t) tree_pos.out.stop; k++) {
+                k != (tsk_size_t) tree_pos.out.stop; k++) {
                 e = out_order[k];
                 p = edge_parent[e];
                 c = edge_child[e];
                 tsk_matvec_calculator_remove_edge(self, p, c);
             }
             for (j = (tsk_size_t) tree_pos.in.start; j != (tsk_size_t) tree_pos.in.stop;
-                 j++) {
+                j++) {
                 e = in_order[j];
                 p = edge_parent[e];
                 c = edge_child[e];
diff --git a/extern/README.md b/extern/README.md
index befd09b..10900dd 100644
--- a/extern/README.md
+++ b/extern/README.md
@@ -1,9 +1,9 @@
-# README for extern/tskit
+# `README` for `extern/tskit`
 
-This directory holds git submodule for tskit and instructions on copying its
-C code into R package.
+This directory holds git submodule for `tskit` and instructions on copying its
+`C` code into `RcppTskit` package.
 
-All commands are run from the root of the RcppTskit repository:
+All commands are run from the root of the `RcppTskit` repository:
 
 ```
 ls -l
@@ -11,9 +11,9 @@ ls -l
 # RcppTskit # R package
 ```
 
-Below instructions show how to add, update, and inspect the tskit code changes.
+Below instructions show how to add, update, and inspect the `tskit` code changes.
 
-## tskit
+## `tskit`
 
 Set up git submodule for the first time:
 
@@ -51,9 +51,9 @@ cat extern/tskit/LICENSE
 # MIT License
 ls -l extern/tskit/c
 cat extern/tskit/c/VERSION.txt
-# 1.3.0
+# 1.3.1
 less extern/tskit/c/CHANGELOG.rst
-# [1.3.0] - 2025-11-27
+# [1.3.1] - 2026-03-06
 ```
 
 Commit git submodule change:
@@ -61,7 +61,7 @@ Commit git submodule change:
 ```
 git status
 git add extern/tskit
-git commit -m "Update tskit submodule to C API version 1.3.0"
+git commit -m "Update tskit submodule to C API version 1.3.1"
 git push
 ```
 
@@ -163,9 +163,10 @@ cp -i extern/tskit/c/subprojects/kastore/kastore.c RcppTskit/src/tskit/.
 Commit changes:
 
 ```
+git diff
 git status
 git add RcppTskit/inst/include/tskit
 git add RcppTskit/src
-git commit -m "Update tskit C API to version 1.3.0"
+git commit -m "Update tskit C API to version 1.3.1"
 git push
 ```
diff --git a/extern/tskit b/extern/tskit
index afacd4f..62f24f4 160000
--- a/extern/tskit
+++ b/extern/tskit
@@ -1 +1 @@
-Subproject commit afacd4f21412c22dc036bc53b97d0ebdab5ae98b
+Subproject commit 62f24f4ca7d1ec91e7c8095aba310a0d7d78b191