diff --git a/.travis.yml b/.travis.yml
index d9936c96f8..3565a3c748 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -8,8 +8,6 @@ language: python
 
 matrix:
   include:
-    - os: linux
-      python: 2.7
     - os: linux
       python: 3.6
     - os: osx
@@ -21,8 +19,6 @@ install:
   - if [[ "$TRAVIS_OS_NAME" == "osx" ]]; then
       curl https://repo.continuum.io/miniconda/Miniconda3-latest-MacOSX-x86_64.sh > miniconda.sh;
       export TRAVIS_PYTHON_VERSION="3.6";
-    elif [[ "$TRAVIS_PYTHON_VERSION" == "2.7" ]]; then
-      wget https://repo.continuum.io/miniconda/Miniconda2-latest-Linux-x86_64.sh -O miniconda.sh;
     else
       wget https://repo.continuum.io/miniconda/Miniconda3-latest-Linux-x86_64.sh -O miniconda.sh;
     fi
diff --git a/docs/Makefile b/docs/Makefile
index a1d061f2df..3043878fc1 100644
--- a/docs/Makefile
+++ b/docs/Makefile
@@ -3,7 +3,7 @@
 
 # You can set these variables from the command line.
 SPHINXOPTS    = -W
-SPHINXBUILD   = sphinx-build
+SPHINXBUILD   = python3 -m sphinx
 SOURCEDIR     = .
 BUILDDIR      = _build
 
diff --git a/docs/conf.py b/docs/conf.py
index b68ac878be..ae734a5e41 100644
--- a/docs/conf.py
+++ b/docs/conf.py
@@ -43,7 +43,7 @@ def __getattr__(cls, name):
 # -- Project information -----------------------------------------------------
 
 project = 'tskit'
-copyright = '2018, Tskit developers'
+copyright = '2018-2019, Tskit developers'
 author = 'Tskit developers'
 
 # The short X.Y version
@@ -76,7 +76,7 @@ def handle_item(fieldarg, content):
             # inconsistencies later when references are resolved
             fieldtype = types.pop(fieldarg)
             if len(fieldtype) == 1 and isinstance(fieldtype[0], nodes.Text):
-                typename = u''.join(n.astext() for n in fieldtype)
+                typename = ''.join(n.astext() for n in fieldtype)
                 par.extend(self.make_xrefs(self.typerolename, domain, typename,
                                            addnodes.literal_emphasis))
             else:
@@ -117,8 +117,12 @@ def handle_item(fieldarg, content):
     'sphinx.ext.intersphinx',
     'sphinx.ext.todo',
     'breathe',
+    'sphinx_issues',
 ]
 
+# Github repo
+issues_github_path = "tskit-dev/tskit"
+
 # Add any paths that contain templates here, relative to this directory.
 templates_path = ['_templates']
 
diff --git a/docs/development.rst b/docs/development.rst
index ad5f502c7b..6445c4baf3 100644
--- a/docs/development.rst
+++ b/docs/development.rst
@@ -25,8 +25,7 @@ Quickstart
 
 - Install the Python development requirements using
   ``pip install -r python/requirements/development.txt``.
-- Build the low level module by running ``make`` in the ``python`` directory. Python 3.x
-  is the default for developement (Python 2.x is discouraged).
+- Build the low level module by running ``make`` in the ``python`` directory.
 - Run the tests to ensure everything has worked: ``python -m nose -vs``. These should
   all pass.
 - Make your changes in a local branch, and open a pull request on GitHub when you
diff --git a/docs/rtd_requirements.txt b/docs/rtd_requirements.txt
index 4f9895802f..667bb2d02a 100644
--- a/docs/rtd_requirements.txt
+++ b/docs/rtd_requirements.txt
@@ -1,6 +1,6 @@
 numpy
-six
 svgwrite
 jsonschema
 h5py
 breathe
+sphinx-issues
diff --git a/python/CHANGELOG.rst b/python/CHANGELOG.rst
index e9797a61c0..5fabe3960d 100644
--- a/python/CHANGELOG.rst
+++ b/python/CHANGELOG.rst
@@ -1,3 +1,15 @@
+--------------------
+[0.1.5] - 2019-XX-XX
+--------------------
+
+Under development.
+
+This release removes support for Python 2.
+
+**New features**
+
+- Remove support for Python 2 (:user:`hugovk`). See :issue:`137` and :pr:`140`.
+
 --------------------
 [0.1.4] - 2019-02-01
 --------------------
diff --git a/python/Makefile b/python/Makefile
index 3b9e7e3f5f..32ec4b43a8 100644
--- a/python/Makefile
+++ b/python/Makefile
@@ -8,9 +8,6 @@ allchecks: _tskitmodule.c
 ext3: _tskitmodule.c 
 	python3 setup.py build_ext --inplace
 
-ext2: _tskitmodule.c 
-	python2 setup.py build_ext --inplace
-
 ctags:
 	ctags lib/*.c lib/*.h tskit/*.py
 
diff --git a/python/_tskitmodule.c b/python/_tskitmodule.c
index b0bf126e1f..d40e98fc6d 100644
--- a/python/_tskitmodule.c
+++ b/python/_tskitmodule.c
@@ -1,21 +1,27 @@
 /*
-** Copyright (C) 2014-2018 University of Oxford
-**
-** This file is part of tskit.
-**
-** tskit is free software: you can redistribute it and/or modify
-** it under the terms of the GNU General Public License as published by
-** the Free Software Foundation, either version 3 of the License, or
-** (at your option) any later version.
-**
-** tskit is distributed in the hope that it will be useful,
-** but WITHOUT ANY WARRANTY; without even the implied warranty of
-** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
-** GNU General Public License for more details.
-**
-** You should have received a copy of the GNU General Public License
-** along with tskit.  If not, see <http://www.gnu.org/licenses/>.
-*/
+ * MIT License
+ *
+ * Copyright (c) 2019 Tskit Developers
+ * Copyright (c) 2015-2018 University of Oxford
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a copy
+ * of this software and associated documentation files (the "Software"), to deal
+ * in the Software without restriction, including without limitation the rights
+ * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+ * copies of the Software, and to permit persons to whom the Software is
+ * furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included in all
+ * copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+ * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+ * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+ * SOFTWARE.
+ */
 
 #define PY_SSIZE_T_CLEAN
 #define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION
@@ -28,13 +34,6 @@
 #include "kastore.h"
 #include "tskit.h"
 
-#if PY_MAJOR_VERSION >= 3
-#define IS_PY3K
-#endif
-
-#define MODULE_DOC \
-"Low level interface for tskit"
-
 #define SET_COLS 0
 #define APPEND_COLS 1
 
@@ -8427,8 +8426,6 @@ static PyTypeObject LdCalculatorType = {
 static PyObject *
 tskit_get_kastore_version(PyObject *self)
 {
-    /* TODO if we provide the option of linking against kastore separately, we
-     * should return the link time version using kas_get_version */
     return Py_BuildValue("iii", KAS_VERSION_MAJOR, KAS_VERSION_MINOR, KAS_VERSION_PATCH);
 }
 
@@ -8446,49 +8443,28 @@ static PyMethodDef tskit_methods[] = {
     {NULL}        /* Sentinel */
 };
 
-/* Initialisation code supports Python 2.x and 3.x. The framework uses the
- * recommended structure from http://docs.python.org/howto/cporting.html.
- * I've ignored the point about storing state in globals, as the examples
- * from the Python documentation still use this idiom.
- */
-
-#if PY_MAJOR_VERSION >= 3
-
 static struct PyModuleDef tskitmodule = {
     PyModuleDef_HEAD_INIT,
-    "_tskit",   /* name of module */
-    MODULE_DOC, /* module documentation, may be NULL */
+    "_tskit",
+    "Low level interface for tskit",
     -1,
     tskit_methods,
     NULL, NULL, NULL, NULL
 };
 
-#define INITERROR return NULL
-
 PyObject *
 PyInit__tskit(void)
-
-#else
-#define INITERROR return
-
-void
-init_tskit(void)
-#endif
 {
-#if PY_MAJOR_VERSION >= 3
     PyObject *module = PyModule_Create(&tskitmodule);
-#else
-    PyObject *module = Py_InitModule3("_tskit", tskit_methods, MODULE_DOC);
-#endif
     if (module == NULL) {
-        INITERROR;
+        return NULL;
     }
     import_array();
 
     /* LightweightTableCollection type */
     LightweightTableCollectionType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&LightweightTableCollectionType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&LightweightTableCollectionType);
     PyModule_AddObject(module, "LightweightTableCollection",
@@ -8497,7 +8473,7 @@ init_tskit(void)
     /* IndividualTable type */
     IndividualTableType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&IndividualTableType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&IndividualTableType);
     PyModule_AddObject(module, "IndividualTable", (PyObject *) &IndividualTableType);
@@ -8505,7 +8481,7 @@ init_tskit(void)
     /* NodeTable type */
     NodeTableType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&NodeTableType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&NodeTableType);
     PyModule_AddObject(module, "NodeTable", (PyObject *) &NodeTableType);
@@ -8513,7 +8489,7 @@ init_tskit(void)
     /* EdgeTable type */
     EdgeTableType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&EdgeTableType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&EdgeTableType);
     PyModule_AddObject(module, "EdgeTable", (PyObject *) &EdgeTableType);
@@ -8521,7 +8497,7 @@ init_tskit(void)
     /* MigrationTable type */
     MigrationTableType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&MigrationTableType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&MigrationTableType);
     PyModule_AddObject(module, "MigrationTable", (PyObject *) &MigrationTableType);
@@ -8529,7 +8505,7 @@ init_tskit(void)
     /* SiteTable type */
     SiteTableType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&SiteTableType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&SiteTableType);
     PyModule_AddObject(module, "SiteTable", (PyObject *) &SiteTableType);
@@ -8537,7 +8513,7 @@ init_tskit(void)
     /* MutationTable type */
     MutationTableType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&MutationTableType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&MutationTableType);
     PyModule_AddObject(module, "MutationTable", (PyObject *) &MutationTableType);
@@ -8545,7 +8521,7 @@ init_tskit(void)
     /* PopulationTable type */
     PopulationTableType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&PopulationTableType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&PopulationTableType);
     PyModule_AddObject(module, "PopulationTable", (PyObject *) &PopulationTableType);
@@ -8553,7 +8529,7 @@ init_tskit(void)
     /* ProvenanceTable type */
     ProvenanceTableType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&ProvenanceTableType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&ProvenanceTableType);
     PyModule_AddObject(module, "ProvenanceTable", (PyObject *) &ProvenanceTableType);
@@ -8561,7 +8537,7 @@ init_tskit(void)
     /* TableCollectionTable type */
     TableCollectionType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&TableCollectionType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&TableCollectionType);
     PyModule_AddObject(module, "TableCollection", (PyObject *) &TableCollectionType);
@@ -8569,7 +8545,7 @@ init_tskit(void)
     /* TreeSequence type */
     TreeSequenceType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&TreeSequenceType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&TreeSequenceType);
     PyModule_AddObject(module, "TreeSequence", (PyObject *) &TreeSequenceType);
@@ -8577,7 +8553,7 @@ init_tskit(void)
     /* Tree type */
     TreeType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&TreeType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&TreeType);
     PyModule_AddObject(module, "Tree", (PyObject *) &TreeType);
@@ -8585,7 +8561,7 @@ init_tskit(void)
     /* TreeDiffIterator type */
     TreeDiffIteratorType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&TreeDiffIteratorType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&TreeDiffIteratorType);
     PyModule_AddObject(module, "TreeDiffIterator", (PyObject *) &TreeDiffIteratorType);
@@ -8593,7 +8569,7 @@ init_tskit(void)
     /* VcfConverter type */
     VcfConverterType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&VcfConverterType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&VcfConverterType);
     PyModule_AddObject(module, "VcfConverter", (PyObject *) &VcfConverterType);
@@ -8601,7 +8577,7 @@ init_tskit(void)
     /* HaplotypeGenerator type */
     HaplotypeGeneratorType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&HaplotypeGeneratorType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&HaplotypeGeneratorType);
     PyModule_AddObject(module, "HaplotypeGenerator",
@@ -8610,7 +8586,7 @@ init_tskit(void)
     /* VariantGenerator type */
     VariantGeneratorType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&VariantGeneratorType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&VariantGeneratorType);
     PyModule_AddObject(module, "VariantGenerator", (PyObject *) &VariantGeneratorType);
@@ -8618,7 +8594,7 @@ init_tskit(void)
     /* LdCalculator type */
     LdCalculatorType.tp_new = PyType_GenericNew;
     if (PyType_Ready(&LdCalculatorType) < 0) {
-        INITERROR;
+        return NULL;
     }
     Py_INCREF(&LdCalculatorType);
     PyModule_AddObject(module, "LdCalculator", (PyObject *) &LdCalculatorType);
@@ -8652,7 +8628,5 @@ init_tskit(void)
     PyModule_AddIntConstant(module, "FORWARD", TSK_DIR_FORWARD);
     PyModule_AddIntConstant(module, "REVERSE", TSK_DIR_REVERSE);
 
-#if PY_MAJOR_VERSION >= 3
     return module;
-#endif
 }
diff --git a/python/requirements/conda-minimal.txt b/python/requirements/conda-minimal.txt
index 9ed290c21e..6f4570efc1 100644
--- a/python/requirements/conda-minimal.txt
+++ b/python/requirements/conda-minimal.txt
@@ -1,7 +1,6 @@
 numpy
 nose
 h5py
-six
 jsonschema
 svgwrite
 msprime
diff --git a/python/requirements/development.txt b/python/requirements/development.txt
index 6674bd2c25..2f62b5fd9f 100644
--- a/python/requirements/development.txt
+++ b/python/requirements/development.txt
@@ -6,12 +6,12 @@ mock
 newick
 nose
 numpy
-six
 kastore
 jsonschema
 sphinx
 sphinx-argparse
 sphinx_rtd_theme
+sphinx-issues
 breathe
 svgwrite
 pyparsing 
diff --git a/python/setup.py b/python/setup.py
index 9db138ff71..ccbdbc13a7 100644
--- a/python/setup.py
+++ b/python/setup.py
@@ -1,4 +1,3 @@
-import sys
 import os.path
 import codecs
 import platform
@@ -14,10 +13,7 @@
 class local_build_ext(build_ext):
     def finalize_options(self):
         build_ext.finalize_options(self)
-        if sys.version_info[0] >= 3:
-            import builtins
-        else:
-            import __builtin__ as builtins
+        import builtins
         # Prevent numpy from thinking it is still in its setup process:
         builtins.__NUMPY_SETUP__ = False
         import numpy
@@ -77,23 +73,23 @@ def finalize_options(self):
     version=tskit_version,
     # TODO setup a tskit developers email address.
     author_email='jerome.kelleher@well.ox.ac.uk',
+    python_requires='>=3.5',
     classifiers=[
         'Development Status :: 4 - Beta',
         'Intended Audience :: Developers',
         'Topic :: Scientific/Engineering :: Bio-Informatics',
         'License :: OSI Approved :: MIT License',
-        'Programming Language :: Python :: 2',
-        'Programming Language :: Python :: 2.7',
         'Programming Language :: Python :: 3',
         'Programming Language :: Python :: 3.5',
         'Programming Language :: Python :: 3.6',
         'Programming Language :: Python :: 3.7',
+        'Programming Language :: Python :: 3 :: Only',
     ],
     keywords='tree sequence',
     packages=['tskit'],
     include_package_data=True,
     ext_modules=[_tskit_module],
-    install_requires=[numpy_ver, "h5py", "svgwrite", "six", "jsonschema"],
+    install_requires=[numpy_ver, "h5py", "svgwrite", "jsonschema"],
     entry_points={
         'console_scripts': [
             'tskit=tskit.__main__:main',
diff --git a/python/stress_lowlevel.py b/python/stress_lowlevel.py
index ce5b25d40b..0f931c90c6 100644
--- a/python/stress_lowlevel.py
+++ b/python/stress_lowlevel.py
@@ -2,9 +2,6 @@
 Code to stress the low-level API as much as possible to expose
 any memory leaks or error handling issues.
 """
-from __future__ import print_function
-from __future__ import division
-
 import argparse
 import unittest
 import random
diff --git a/python/tests/__init__.py b/python/tests/__init__.py
index ca66d87f3c..6bced163ba 100644
--- a/python/tests/__init__.py
+++ b/python/tests/__init__.py
@@ -20,10 +20,6 @@
 # OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
 # SOFTWARE.
 
-from __future__ import print_function
-from __future__ import unicode_literals
-from __future__ import division
-
 import base64
 
 # TODO remove this code and refactor elsewhere.
diff --git a/python/tests/simplify.py b/python/tests/simplify.py
index 2d5759258c..0c8df77a85 100644
--- a/python/tests/simplify.py
+++ b/python/tests/simplify.py
@@ -23,10 +23,6 @@
 """
 Python implementation of the simplify algorithm.
 """
-from __future__ import print_function
-from __future__ import unicode_literals
-from __future__ import division
-
 import sys
 
 import numpy as np
@@ -288,7 +284,7 @@ def process_parent_edges(self, edges):
         """
         Process all of the edges for a given parent.
         """
-        assert len(set(e.parent for e in edges)) == 1
+        assert len({e.parent for e in edges}) == 1
         parent = edges[0].parent
         S = []
         for edge in edges:
diff --git a/python/tests/test_dict_encoding.py b/python/tests/test_dict_encoding.py
index f3a7a10d62..774c81834b 100644
--- a/python/tests/test_dict_encoding.py
+++ b/python/tests/test_dict_encoding.py
@@ -23,9 +23,6 @@
 Test cases for the low-level dictionary encoding used to move
 data around in C.
 """
-from __future__ import print_function
-from __future__ import division
-
 import unittest
 
 import msprime
@@ -138,7 +135,7 @@ def test_missing_sequence_length(self):
     def test_missing_tables(self):
         tables = get_example_tables()
         d = tables.asdict()
-        table_names = set(d.keys()) - set(["sequence_length"])
+        table_names = set(d.keys()) - {"sequence_length"}
         for table_name in table_names:
             d = tables.asdict()
             del d[table_name]
@@ -149,7 +146,7 @@ def test_missing_tables(self):
     def test_missing_columns(self):
         tables = get_example_tables()
         d = tables.asdict()
-        table_names = set(d.keys()) - set(["sequence_length"])
+        table_names = set(d.keys()) - {"sequence_length"}
         for table_name in table_names:
             table_dict = d[table_name]
             for colname in table_dict.keys():
@@ -170,7 +167,7 @@ class TestBadTypes(unittest.TestCase):
     def verify_columns(self, value):
         tables = get_example_tables()
         d = tables.asdict()
-        table_names = set(d.keys()) - set(["sequence_length"])
+        table_names = set(d.keys()) - {"sequence_length"}
         for table_name in table_names:
             table_dict = d[table_name]
             for colname in table_dict.keys():
@@ -208,7 +205,7 @@ def verify(self, num_rows):
 
         tables = get_example_tables()
         d = tables.asdict()
-        table_names = set(d.keys()) - set(["sequence_length"])
+        table_names = set(d.keys()) - {"sequence_length"}
         for table_name in sorted(table_names):
             table_dict = d[table_name]
             for colname in sorted(table_dict.keys()):
diff --git a/python/tests/test_drawing.py b/python/tests/test_drawing.py
index 5b6e13a5ab..c8099b6ff7 100644
--- a/python/tests/test_drawing.py
+++ b/python/tests/test_drawing.py
@@ -24,22 +24,16 @@
 """
 Test cases for visualisation in tskit.
 """
-from __future__ import print_function
-from __future__ import division
-
+import io
 import os
-import sys
 import tempfile
 import unittest
 import xml.etree
 
 import msprime
-import six
 import tskit
 import tests.tsutil as tsutil
 
-IS_PY2 = sys.version_info[0] < 3
-
 
 class TestTreeDraw(unittest.TestCase):
     """
@@ -169,13 +163,10 @@ def test_ascii_variants(self):
     def test_unicode_variants(self):
         t = self.get_binary_tree()
         for fmt in ["unicode", "UNICODE", "uniCODE"]:
-            if IS_PY2:
-                self.assertRaises(ValueError, t.draw, format=fmt)
-            else:
-                output = t.draw(format=fmt)
-                self.assertRaises(
-                    xml.etree.ElementTree.ParseError, xml.etree.ElementTree.fromstring,
-                    output)
+            output = t.draw(format=fmt)
+            self.assertRaises(
+                xml.etree.ElementTree.ParseError, xml.etree.ElementTree.fromstring,
+                output)
 
     def test_bad_formats(self):
         t = self.get_binary_tree()
@@ -236,7 +227,7 @@ def test_draw_zero_edge_tree(self):
         self.verify_basic_text(text)
 
     def test_even_num_children_tree(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           1
@@ -246,7 +237,7 @@ def test_even_num_children_tree(self):
         5   1           5
         6   1           7
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       6       0
         0       1       6       1
@@ -261,7 +252,7 @@ def test_even_num_children_tree(self):
         self.verify_basic_text(text)
 
     def test_odd_num_children_tree(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           1
@@ -270,7 +261,7 @@ def test_odd_num_children_tree(self):
         4   1           4
         5   1           5
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       5       0
         0       1       5       1
@@ -302,7 +293,6 @@ def test_no_node_labels(self):
             self.assertEqual(text.find(str(u)), -1)
 
 
-@unittest.skipIf(IS_PY2, "Unicode tree drawing not supported on Python 2")
 class TestDrawUnicode(TestDrawText):
     """
     Tests the Unicode tree drawing method
@@ -324,13 +314,13 @@ def verify_text_rendering(self, drawn, drawn_tree, debug=False):
             self.assertEqual(l1.rstrip(), l2.rstrip())
 
     def test_simple_tree(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           0
         2   1           2
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       2       0
         0       1       2       1
@@ -345,14 +335,14 @@ def test_simple_tree(self):
         self.verify_text_rendering(drawn, tree)
 
     def test_trident_tree(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           0
         2   1           0
         3   1           2
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       3       0
         0       1       3       1
@@ -368,7 +358,7 @@ def test_trident_tree(self):
         self.verify_text_rendering(drawn, tree)
 
     def test_pitchfork_tree(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           0
@@ -376,7 +366,7 @@ def test_pitchfork_tree(self):
         3   1           0
         4   1           2
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       4       0
         0       1       4       1
@@ -405,13 +395,13 @@ def test_pitchfork_tree(self):
         self.verify_text_rendering(drawn, tree)
 
     def test_stick_tree(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           1
         2   1           2
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       1       0
         0       1       2       1
diff --git a/python/tests/test_file_format.py b/python/tests/test_file_format.py
index 018de57450..cd9775e436 100644
--- a/python/tests/test_file_format.py
+++ b/python/tests/test_file_format.py
@@ -23,9 +23,6 @@
 """
 Test cases for tskit's file format.
 """
-from __future__ import print_function
-from __future__ import division
-
 import os
 import tempfile
 import unittest
diff --git a/python/tests/test_highlevel.py b/python/tests/test_highlevel.py
index 2d96a661ad..fd6d41b20e 100644
--- a/python/tests/test_highlevel.py
+++ b/python/tests/test_highlevel.py
@@ -23,29 +23,18 @@
 """
 Test cases for the high level interface to msprime.
 """
-from __future__ import print_function
-from __future__ import division
-
-try:
-    # We use the zip as iterator functionality here.
-    from future_builtins import zip
-except ImportError:
-    # This fails for Python 3.x, but that's fine.
-    pass
-
 import collections
 import itertools
+import io
 import json
 import math
 import os
 import random
 import shutil
-import six
 import tempfile
 import unittest
 import warnings
 import uuid as _uuid
-import sys
 
 import numpy as np
 import msprime
@@ -56,8 +45,6 @@
 import tests.tsutil as tsutil
 import tests.simplify as simplify
 
-IS_PY2 = sys.version_info[0] < 3
-
 
 def insert_uniform_mutations(tables, num_mutations, nodes):
     """
@@ -1105,7 +1092,7 @@ def verify_samples(self, ts):
     def test_samples(self):
         for ts in get_example_tree_sequences():
             self.verify_samples(ts)
-            pops = set(node.population for node in ts.nodes())
+            pops = {node.population for node in ts.nodes()}
             for pop in pops:
                 subsample = ts.samples(pop)
                 self.assertTrue(np.array_equal(subsample, ts.samples(population=pop)))
@@ -1361,7 +1348,7 @@ def test_simplify(self):
             num_mutations += ts.get_num_mutations()
             sample_sizes = {0, 1}
             if n > 2:
-                sample_sizes |= set([2, max(2, n // 2), n - 1])
+                sample_sizes |= {2, max(2, n // 2), n - 1}
             for k in sample_sizes:
                 subset = random.sample(list(ts.samples()), k)
                 self.verify_simplify_topology(ts, subset)
@@ -1718,10 +1705,10 @@ def convert(v):
     def test_output_format(self):
         for ts in get_example_tree_sequences():
             for precision in [2, 7]:
-                nodes_file = six.StringIO()
-                edges_file = six.StringIO()
-                sites_file = six.StringIO()
-                mutations_file = six.StringIO()
+                nodes_file = io.StringIO()
+                edges_file = io.StringIO()
+                sites_file = io.StringIO()
+                mutations_file = io.StringIO()
                 ts.dump_text(
                     nodes=nodes_file, edges=edges_file, sites=sites_file,
                     mutations=mutations_file, precision=precision)
@@ -1781,12 +1768,12 @@ def verify_approximate_equality(self, ts1, ts2):
 
     def test_text_record_round_trip(self):
         for ts1 in get_example_tree_sequences():
-            nodes_file = six.StringIO()
-            edges_file = six.StringIO()
-            sites_file = six.StringIO()
-            mutations_file = six.StringIO()
-            individuals_file = six.StringIO()
-            populations_file = six.StringIO()
+            nodes_file = io.StringIO()
+            edges_file = io.StringIO()
+            sites_file = io.StringIO()
+            mutations_file = io.StringIO()
+            individuals_file = io.StringIO()
+            populations_file = io.StringIO()
             ts1.dump_text(
                 nodes=nodes_file, edges=edges_file, sites=sites_file,
                 mutations=mutations_file, individuals=individuals_file,
@@ -1806,20 +1793,20 @@ def test_text_record_round_trip(self):
             self.verify_approximate_equality(ts1, ts2)
 
     def test_empty_files(self):
-        nodes_file = six.StringIO("is_sample\ttime\n")
-        edges_file = six.StringIO("left\tright\tparent\tchild\n")
-        sites_file = six.StringIO("position\tancestral_state\n")
-        mutations_file = six.StringIO("site\tnode\tderived_state\n")
+        nodes_file = io.StringIO("is_sample\ttime\n")
+        edges_file = io.StringIO("left\tright\tparent\tchild\n")
+        sites_file = io.StringIO("position\tancestral_state\n")
+        mutations_file = io.StringIO("site\tnode\tderived_state\n")
         self.assertRaises(
             _tskit.LibraryError, tskit.load_text,
             nodes=nodes_file, edges=edges_file, sites=sites_file,
             mutations=mutations_file)
 
     def test_empty_files_sequence_length(self):
-        nodes_file = six.StringIO("is_sample\ttime\n")
-        edges_file = six.StringIO("left\tright\tparent\tchild\n")
-        sites_file = six.StringIO("position\tancestral_state\n")
-        mutations_file = six.StringIO("site\tnode\tderived_state\n")
+        nodes_file = io.StringIO("is_sample\ttime\n")
+        edges_file = io.StringIO("left\tright\tparent\tchild\n")
+        sites_file = io.StringIO("position\tancestral_state\n")
+        mutations_file = io.StringIO("site\tnode\tderived_state\n")
         ts = tskit.load_text(
             nodes=nodes_file, edges=edges_file, sites=sites_file,
             mutations=mutations_file, sequence_length=100)
@@ -2336,8 +2323,8 @@ def verify_random_permutation(self, ts):
         other_ts.dump(self.temp_file)
         ts3 = tskit.load(self.temp_file)
         self.verify_tree_sequences_equal(other_ts, ts3)
-        nodes_file = six.StringIO()
-        edges_file = six.StringIO()
+        nodes_file = io.StringIO()
+        edges_file = io.StringIO()
         # Also verify we can read the text version.
         other_ts.dump_text(nodes=nodes_file, edges=edges_file, precision=14)
         nodes_file.seek(0)
@@ -2414,7 +2401,6 @@ def get_instances(self, n):
             tskit.Edge(left=j, right=j, parent=j, child=j) for j in range(n)]
 
 
-@unittest.skipIf("py2", IS_PY2)
 class TestSiteContainer(unittest.TestCase, SimpleContainersMixin):
     def get_instances(self, n):
         return [
@@ -2458,7 +2444,6 @@ def get_instances(self, n):
             tskit.Edgeset(left=j, right=j, parent=j, children=j) for j in range(n)]
 
 
-@unittest.skipIf("py2", IS_PY2)
 class TestVariantContainer(unittest.TestCase, SimpleContainersMixin):
     def get_instances(self, n):
         return [
diff --git a/python/tests/test_lowlevel.py b/python/tests/test_lowlevel.py
index 49ab96bb61..64276483a0 100644
--- a/python/tests/test_lowlevel.py
+++ b/python/tests/test_lowlevel.py
@@ -23,16 +23,11 @@
 """
 Test cases for the low level C interface to tskit.
 """
-from __future__ import print_function
-from __future__ import division
-from __future__ import unicode_literals
-
 import collections
 import itertools
 import os
 import platform
 import random
-import sys
 import tempfile
 import unittest
 
@@ -40,7 +35,6 @@
 
 import _tskit
 
-IS_PY2 = sys.version_info[0] < 3
 IS_WINDOWS = platform.system() == "Windows"
 
 
diff --git a/python/tests/test_metadata.py b/python/tests/test_metadata.py
index 69c1c1b820..f72fa0af1e 100644
--- a/python/tests/test_metadata.py
+++ b/python/tests/test_metadata.py
@@ -24,9 +24,7 @@
 """
 Tests for metadata handling.
 """
-from __future__ import print_function
-from __future__ import division
-
+import io
 import json
 import os
 import tempfile
@@ -35,7 +33,6 @@
 
 import numpy as np
 import python_jsonschema_objects as pjs
-import six
 import msprime
 
 import tskit
@@ -192,7 +189,7 @@ class TestLoadTextMetadata(unittest.TestCase):
     """
 
     def test_individuals(self):
-        individuals = six.StringIO("""\
+        individuals = io.StringIO("""\
         id  flags location     metadata
         0   1     0.0,1.0,0.0  abc
         1   1     1.0,2.0      XYZ+
@@ -212,7 +209,7 @@ def test_individuals(self):
                              b.metadata)
 
     def test_nodes(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time    metadata
         0   1           0   abc
         1   1           0   XYZ+
@@ -226,7 +223,7 @@ def test_nodes(self):
                              b.metadata)
 
     def test_sites(self):
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         position    ancestral_state metadata
         0.1 A   abc
         0.5 C   XYZ+
@@ -240,7 +237,7 @@ def test_sites(self):
                              b.metadata)
 
     def test_mutations(self):
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state   metadata
         0   2   C   mno
         0   3   G   )(*&^%$#@!
@@ -253,7 +250,7 @@ def test_mutations(self):
                              b.metadata)
 
     def test_populations(self):
-        populations = six.StringIO("""\
+        populations = io.StringIO("""\
         id    metadata
         0     mno
         1     )(*&^%$#@!
diff --git a/python/tests/test_newick.py b/python/tests/test_newick.py
index e5d3e47487..a6bbfc8c5c 100644
--- a/python/tests/test_newick.py
+++ b/python/tests/test_newick.py
@@ -23,9 +23,6 @@
 """
 Tests for the newick output feature.
 """
-from __future__ import print_function
-from __future__ import division
-
 import unittest
 
 import msprime
diff --git a/python/tests/test_provenance.py b/python/tests/test_provenance.py
index faaf0078d8..3da9d675d8 100644
--- a/python/tests/test_provenance.py
+++ b/python/tests/test_provenance.py
@@ -23,10 +23,6 @@
 """
 Tests for the provenance information attached to tree sequences.
 """
-from __future__ import print_function
-from __future__ import unicode_literals
-from __future__ import division
-
 import unittest
 import json
 import platform
diff --git a/python/tests/test_stats.py b/python/tests/test_stats.py
index 709ef6f917..8a33d4db9d 100644
--- a/python/tests/test_stats.py
+++ b/python/tests/test_stats.py
@@ -23,11 +23,7 @@
 """
 Test cases for stats calculations in tskit.
 """
-from __future__ import print_function
-from __future__ import division
-
 import unittest
-import sys
 
 import numpy as np
 import msprime
@@ -38,9 +34,6 @@
 import tests.test_wright_fisher as wf
 
 
-IS_PY2 = sys.version_info[0] < 3
-
-
 def get_r2_matrix(ts):
     """
     Returns the matrix for the specified tree sequence. This is computed
@@ -377,14 +370,8 @@ def verify(self, ts, reference_sets, focal=None):
         A1 = naive_genealogical_nearest_neighbours(ts, focal, reference_sets)
         A2 = tsutil.genealogical_nearest_neighbours(ts, focal, reference_sets)
         A3 = ts.genealogical_nearest_neighbours(focal, reference_sets)
-        if IS_PY2:
-            # Threads not supported on PY2
-            self.assertRaises(
-                ValueError, ts.genealogical_nearest_neighbours, focal,
-                reference_sets, num_threads=3)
-        else:
-            A4 = ts.genealogical_nearest_neighbours(focal, reference_sets, num_threads=3)
-            self.assertTrue(np.array_equal(A3, A4))
+        A4 = ts.genealogical_nearest_neighbours(focal, reference_sets, num_threads=3)
+        self.assertTrue(np.array_equal(A3, A4))
         self.assertEqual(A1.shape, A2.shape)
         self.assertEqual(A1.shape, A3.shape)
         self.assertTrue(np.allclose(A1, A2))
diff --git a/python/tests/test_tables.py b/python/tests/test_tables.py
index 17f76a939e..28b83ff5d8 100644
--- a/python/tests/test_tables.py
+++ b/python/tests/test_tables.py
@@ -1,4 +1,3 @@
-# -*- coding: utf-8 -*-
 # MIT License
 #
 # Copyright (c) 2018-2019 Tskit Developers
@@ -25,18 +24,14 @@
 Test cases for the low-level tables used to transfer information
 between simulations and the tree sequence.
 """
-from __future__ import print_function
-from __future__ import division
-
+import io
 import pickle
 import random
 import string
 import unittest
 import warnings
-import sys
 
 import numpy as np
-import six
 
 import tskit
 import _tskit
@@ -44,8 +39,6 @@
 
 import tests.tsutil as tsutil
 
-IS_PY2 = sys.version_info[0] < 3
-
 
 def random_bytes(max_length):
     """
@@ -889,7 +882,7 @@ def test_random_cases(self):
             self.verify_packing(strings)
 
     def test_unicode(self):
-        self.verify_packing([u'abcdé', u'€'])
+        self.verify_packing(['abcdé', '€'])
 
 
 class TestBytePacking(unittest.TestCase):
@@ -1143,20 +1136,20 @@ class TestSortMutations(unittest.TestCase):
     """
 
     def test_sort_mutations_stability(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         position    ancestral_state
         0.1     0
         0.2     0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state   parent
         1       0       1               -1
         1       1       1               -1
@@ -1176,19 +1169,19 @@ def test_sort_mutations_stability(self):
         self.assertEqual(list(mutations.node), [1, 0, 0, 1])
 
     def test_sort_mutations_remap_parent_id(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         position    ancestral_state
         0.1     0
         0.2     0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state   parent
         1       0       1               -1
         1       0       0               0
@@ -1211,18 +1204,18 @@ def test_sort_mutations_remap_parent_id(self):
         self.assertEqual(list(mutations.parent), [-1, 0, 1, -1, 3, 4])
 
     def test_sort_mutations_bad_parent_id(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         position    ancestral_state
         0.1     0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state   parent
         1       0       1               -2
         """)
@@ -1238,7 +1231,6 @@ class TestSimplifyTables(unittest.TestCase):
     """
     random_seed = 42
 
-    @unittest.skipIf(IS_PY2, "Warnings different in Py2")
     def test_deprecated_zero_mutation_sites(self):
         ts = msprime.simulate(10,  mutation_rate=1, random_seed=self.random_seed)
         tables = ts.dump_tables()
diff --git a/python/tests/test_threads.py b/python/tests/test_threads.py
index d7b8f9779b..46a4f8443f 100644
--- a/python/tests/test_threads.py
+++ b/python/tests/test_threads.py
@@ -23,10 +23,6 @@
 """
 Test cases for threading enabled aspects of the API.
 """
-from __future__ import print_function
-from __future__ import division
-
-import sys
 import threading
 import unittest
 import platform
@@ -37,7 +33,6 @@
 import tskit
 import tests.tsutil as tsutil
 
-IS_PY2 = sys.version_info[0] < 3
 IS_WINDOWS = platform.system() == "Windows"
 
 
@@ -141,10 +136,9 @@ def worker(thread_index, results):
             self.assertEqual(results[j][0], m - j - 1)
 
 
-# @unittest.skipIf(IS_PY2, "Cannot test thread support on Py2.")
-# Temporarily skipping these on windows too. See
+# Temporarily skipping these on Windows. See
 # https://github.com/jeromekelleher/tskit/issues/344
-@unittest.skipIf(IS_PY2 or IS_WINDOWS, "Cannot test thread support on Py2.")
+@unittest.skipIf(IS_WINDOWS, "Cannot test thread support on Windows.")
 class TestTables(unittest.TestCase):
     """
     Tests to ensure that attempts to access tables in threads correctly
diff --git a/python/tests/test_topology.py b/python/tests/test_topology.py
index 57a64f778c..7501b23cb3 100644
--- a/python/tests/test_topology.py
+++ b/python/tests/test_topology.py
@@ -23,21 +23,11 @@
 """
 Test cases for the supported topological variations and operations.
 """
-from __future__ import print_function
-from __future__ import division
-
-try:
-    # We use the zip as iterator functionality here.
-    from future_builtins import zip
-except ImportError:
-    # This fails for Python 3.x, but that's fine.
-    pass
-
+import io
 import unittest
 import itertools
 import random
 
-import six
 import numpy as np
 import msprime
 
@@ -442,13 +432,13 @@ def verify_zero_roots(self, ts):
             self.assertEqual(tree.roots, [])
 
     def test_simple_hole(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           0
         2   0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       2       0
         2       3       2       0
@@ -463,13 +453,13 @@ def test_simple_hole(self):
         self.verify_trees(ts, expected)
 
     def test_simple_hole_zero_roots(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   0           0
         1   0           0
         2   0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       2       0
         2       3       2       0
@@ -485,13 +475,13 @@ def test_simple_hole_zero_roots(self):
         self.verify_zero_roots(ts)
 
     def test_initial_gap(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           0
         2   0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         1       2       2       0,1
         """)
@@ -502,13 +492,13 @@ def test_initial_gap(self):
         self.verify_trees(ts, expected)
 
     def test_initial_gap_zero_roots(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   0           0
         1   0           0
         2   0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         1       2       2       0,1
         """)
@@ -520,13 +510,13 @@ def test_initial_gap_zero_roots(self):
         self.verify_zero_roots(ts)
 
     def test_final_gap(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           0
         2   0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       2       2       0,1
         """)
@@ -537,13 +527,13 @@ def test_final_gap(self):
         self.verify_trees(ts, expected)
 
     def test_final_gap_zero_roots(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   0           0
         1   0           0
         2   0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       2       2       0,1
         """)
@@ -555,13 +545,13 @@ def test_final_gap_zero_roots(self):
         self.verify_zero_roots(ts)
 
     def test_initial_and_final_gap(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           0
         2   0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         1       2       2       0,1
         """)
@@ -573,13 +563,13 @@ def test_initial_and_final_gap(self):
         self.verify_trees(ts, expected)
 
     def test_initial_and_final_gap_zero_roots(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   0           0
         1   0           0
         2   0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         1       2       2       0,1
         """)
@@ -599,7 +589,7 @@ class TestTsinferExamples(TopologyTestCase):
     def test_no_last_tree(self):
         # The last tree was not being generated here because of a bug in
         # the low-level tree generation code.
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   population      time
         0       1       -1              3.00000000000000
         1       1       -1              2.00000000000000
@@ -613,7 +603,7 @@ def test_no_last_tree(self):
         9       1       -1              1.00000000000000
         10      1       -1              1.00000000000000
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         id      left            right           parent  child
         0       62291.41659631  79679.17408763  1       5
         1       62291.41659631  62374.60889677  1       6
@@ -672,12 +662,12 @@ class TestRecordSquashing(TopologyTestCase):
     Tests that we correctly squash adjacent equal records together.
     """
     def test_single_record(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       1       0
         1       2       1       0
@@ -751,7 +741,7 @@ class TestUnaryNodes(TopologyTestCase):
     """
     def test_simple_case(self):
         # Simple case where we have n = 2 and some unary nodes.
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -760,7 +750,7 @@ def test_simple_case(self):
         4       0           2
         5       0           3
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       2       0
         0       1       3       1
@@ -774,8 +764,8 @@ def test_simple_case(self):
             sites += "{} 0\n".format(position)
             mutations += "{} {} 1\n".format(j, j)
         ts = tskit.load_text(
-            nodes=nodes, edges=edges, sites=six.StringIO(sites),
-            mutations=six.StringIO(mutations), strict=False)
+            nodes=nodes, edges=edges, sites=io.StringIO(sites),
+            mutations=io.StringIO(mutations), strict=False)
 
         self.assertEqual(ts.sample_size, 2)
         self.assertEqual(ts.num_nodes, 6)
@@ -815,7 +805,7 @@ def test_ladder_tree(self):
         root = num_unary_nodes + 3
         root_time = num_unary_nodes + 3
         edges += "0    1     {}      1,{}\n".format(root, num_unary_nodes + 2)
-        ts = tskit.load_text(six.StringIO(nodes), six.StringIO(edges), strict=False)
+        ts = tskit.load_text(io.StringIO(nodes), io.StringIO(edges), strict=False)
         t = next(ts.trees())
         self.assertEqual(t.mrca(0, 1), root)
         self.assertEqual(t.tmrca(0, 1), root_time)
@@ -892,7 +882,7 @@ class TestGeneralSamples(TopologyTestCase):
     """
     def test_simple_case(self):
         # Simple case where we have n = 3 and samples starting at n.
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       0           2
         1       0           1
@@ -900,19 +890,19 @@ def test_simple_case(self):
         3       1           0
         4       1           0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       1       2,3
         0       1       0       1,4
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         position    ancestral_state
         0.1     0
         0.2     0
         0.3     0
         0.4     0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state
         0       2       1
         1       3       1
@@ -1060,21 +1050,21 @@ def verify_simplify(
         required output.
         """
         ts = tskit.load_text(
-            nodes=six.StringIO(nodes_before),
-            edges=six.StringIO(edges_before),
-            sites=six.StringIO(sites_before) if sites_before is not None else None,
+            nodes=io.StringIO(nodes_before),
+            edges=io.StringIO(edges_before),
+            sites=io.StringIO(sites_before) if sites_before is not None else None,
             mutations=(
-                six.StringIO(mutations_before)
+                io.StringIO(mutations_before)
                 if mutations_before is not None else None),
             strict=False)
         before = ts.dump_tables()
 
         ts = tskit.load_text(
-            nodes=six.StringIO(nodes_after),
-            edges=six.StringIO(edges_after),
-            sites=six.StringIO(sites_after) if sites_after is not None else None,
+            nodes=io.StringIO(nodes_after),
+            edges=io.StringIO(edges_after),
+            sites=io.StringIO(sites_after) if sites_after is not None else None,
             mutations=(
-                six.StringIO(mutations_after)
+                io.StringIO(mutations_after)
                 if mutations_after is not None else None),
             strict=False,
             sequence_length=before.sequence_length)
@@ -1112,8 +1102,8 @@ def test_unsorted_edges(self):
         1       2       2       0,1
         1       2       3       0,1
         """
-        nodes = tskit.parse_nodes(six.StringIO(nodes_before), strict=False)
-        edges = tskit.parse_edges(six.StringIO(edges_before), strict=False)
+        nodes = tskit.parse_nodes(io.StringIO(nodes_before), strict=False)
+        edges = tskit.parse_edges(io.StringIO(edges_before), strict=False)
         # Cannot use load_text here because it calls sort()
         tables = tskit.TableCollection(sequence_length=2)
         tables.nodes.set_columns(**nodes.asdict())
@@ -1431,7 +1421,7 @@ class TestNonSampleExternalNodes(TopologyTestCase):
     """
     def test_simple_case(self):
         # Simplest case where we have n = 2 and external non-sample nodes.
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -1439,18 +1429,18 @@ def test_simple_case(self):
         3       0           0
         4       0           0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       2       0,1,3,4
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         id  position    ancestral_state
         0   0.1         0
         1   0.2         0
         2   0.3         0
         3   0.4         0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state
         0       0       1
         1       1       1
@@ -1513,20 +1503,20 @@ class TestMultipleRoots(TopologyTestCase):
 
     def test_simplest_degenerate_case(self):
         # Simplest case where we have n = 2 and no edges.
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         id  position    ancestral_state
         0   0.1         0
         1   0.2         0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state
         0       0         1
         1       1         1
@@ -1552,7 +1542,7 @@ def test_simplest_degenerate_case(self):
 
     def test_simplest_non_degenerate_case(self):
         # Simplest case where we have n = 4 and two trees.
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -1561,19 +1551,19 @@ def test_simplest_non_degenerate_case(self):
         4       0           1
         5       0           2
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       4       0,1
         0       1       5       2,3
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         id  position    ancestral_state
         0   0.1         0
         1   0.2         0
         2   0.3         0
         3   0.4         0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state
         0       0       1
         1       1       1
@@ -1610,7 +1600,7 @@ def test_simplest_non_degenerate_case(self):
 
     def test_two_reducable_trees(self):
         # We have n = 4 and two trees, with some unary nodes and non-sample leaves
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -1622,14 +1612,14 @@ def test_two_reducable_trees(self):
         7       0           3
         8       0           0   # Non sample leaf
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1      4         0
         0       1      5         1
         0       1      6         4,5
         0       1      7         2,3,8
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         id  position    ancestral_state
         0   0.1         0
         1   0.2         0
@@ -1637,7 +1627,7 @@ def test_two_reducable_trees(self):
         3   0.4         0
         4   0.5         0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state
         0       0       1
         1       1       1
@@ -1682,7 +1672,7 @@ def test_two_reducable_trees(self):
 
     def test_one_reducable_tree(self):
         # We have n = 4 and two trees. One tree is reducable and the other isn't.
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -1694,7 +1684,7 @@ def test_one_reducable_tree(self):
         7       0           3
         8       0           0   # Non sample leaf
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1      4         0
         0       1      5         1
@@ -1722,7 +1712,7 @@ def test_one_reducable_tree(self):
     # so there might be other problems with it.
     def test_mutations_over_roots(self):
         # Mutations over root nodes should be ok when we have multiple roots.
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -1731,13 +1721,13 @@ def test_mutations_over_roots(self):
         4       0           2
         5       0           2
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       3       0,1
         0       1       4       3
         0       1       5       2
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         id  position    ancestral_state
         0   0.1         0
         1   0.2         0
@@ -1746,7 +1736,7 @@ def test_mutations_over_roots(self):
         4   0.5         0
         5   0.6         0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state
         0       0       1
         1       1       1
@@ -1850,7 +1840,7 @@ def test_partial_non_sample_external_nodes(self):
         #
         #          (0.0, 0.2),                 (0.2, 0.8),         (0.8, 1.0)
 
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -1861,7 +1851,7 @@ def test_partial_non_sample_external_nodes(self):
         6       0           0.7
         7       0           1.0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.0     0.2     4       2,3
         0.2     0.8     4       0,2
@@ -1904,7 +1894,7 @@ def test_partial_non_sample_external_nodes_2(self):
         # 0.0 0     1           2     1   0       2     0   1           2
         #
         #          (0.0, 0.2),         (0.2, 0.8),         (0.8, 1.0)
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -1915,7 +1905,7 @@ def test_partial_non_sample_external_nodes_2(self):
         6       0           1.0
         7       0           0    # Non sample leaf
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.0     0.2     3       2,7
         0.2     0.8     3       0,2
@@ -1961,7 +1951,7 @@ def test_single_offspring_records(self):
         # 0.0    0     1           2           1   0       2          0   1           2
         #
         #          (0.0, 0.2),               (0.2, 0.8),              (0.8, 1.0)
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id  is_sample   time
         0   1           0
         1   1           0
@@ -1972,7 +1962,7 @@ def test_single_offspring_records(self):
         6   0           0.7
         7   0           1.0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.0     0.2     4       2,3
         0.2     0.8     4       0,2
@@ -2056,7 +2046,7 @@ def test_many_single_offspring(self):
             {0: 9, 1: 10, 2: 5, 3: -1, 4: 3, 5: 3, 6: 3, 7: 5, 8: 7, 9: 6, 10: 8}
         ]
         true_haplotypes = ['0100', '0001', '1110']
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -2070,7 +2060,7 @@ def test_many_single_offspring(self):
         9       0           1
         10      0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.5     1.0     10      1
         0.0     0.4     10      2
@@ -2091,14 +2081,14 @@ def test_many_single_offspring(self):
         0.1     0.9     3       4,5
         0.0     0.1     3       4,5,7
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         position    ancestral_state
         0.05        0
         0.15        0
         0.25        0
         0.4         0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state   parent
         0       7       1               -1
         0      10       0               0
@@ -2151,7 +2141,7 @@ def test_tricky_switches(self):
         # 0   1 2   3 4  5  . 0   1 2   3 4  5  . 4   5 2   3 0  1  .
         #                   .                   .                   .
         # 0.0              0.4                 0.5                 1.0
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -2167,7 +2157,7 @@ def test_tricky_switches(self):
         11      0           3
         12      0           4
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left right parent child
         0.0  0.5   6      0
         0.0  0.5   6      1
@@ -2236,7 +2226,7 @@ def test_tricky_simplify(self):
         #  . 0   1 2   3 4  5  . 4   5 2   3 0  1  .
         #  .                   .                   .
         # 0.4                 0.5                 1.0
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -2255,7 +2245,7 @@ def test_tricky_simplify(self):
         14      0           1
         15      0           2
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left right parent child
         0.0  0.5   6      0,1
         0.5  1.0   6      4,5
@@ -2318,7 +2308,7 @@ def test_ancestral_samples(self):
         #              *           *          *           *         *           *
         #          (0.0, 0.2),                 (0.2, 0.8),         (0.8, 1.0)
 
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       0           0
         1       1           0
@@ -2330,7 +2320,7 @@ def test_ancestral_samples(self):
         7       0           1.0
         8       0           0.8
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.0     0.2     4       2,3
         0.2     0.8     4       0,2
@@ -2396,7 +2386,7 @@ def test_all_ancestral_samples(self):
         #
         #          (0.0, 0.2),                 (0.2, 0.8),         (0.8, 1.0)
 
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       0           0
         1       1           0.1
@@ -2408,7 +2398,7 @@ def test_all_ancestral_samples(self):
         7       0           1.0
         8       0           0.8
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.0     0.2     4       2,3
         0.2     0.8     4       0,2
@@ -2454,7 +2444,7 @@ def test_internal_sampled_node(self):
         # 0.0    0      *         *            *  0      *      0    *         *
         #
         #          (0.0, 0.2),                 (0.2, 0.8),         (0.8, 1.0)
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       0           0
         1       1           0.1
@@ -2466,7 +2456,7 @@ def test_internal_sampled_node(self):
         7       0           1.0
         8       0           0.8
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.0     0.2     4       2,3
         0.2     0.8     4       0,2
@@ -2540,14 +2530,14 @@ class TestBadTrees(unittest.TestCase):
     try to create trees.
     """
     def test_simplest_contradictory_children(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
         2       0           1
         3       0           2
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.0     1.0     2       0
         0.0     1.0     3       0
@@ -2556,14 +2546,14 @@ def test_simplest_contradictory_children(self):
         self.assertRaises(_tskit.LibraryError, list, ts.trees())
 
     def test_partial_overlap_contradictory_children(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
         2       0           1
         3       0           2
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.0     1.0     2       0,1
         0.5     1.0     3       0
@@ -2744,8 +2734,8 @@ def test_many_trees(self):
 
     def test_small_tree_internal_samples(self):
         ts = tskit.load_text(
-            nodes=six.StringIO(self.small_tree_ex_nodes),
-            edges=six.StringIO(self.small_tree_ex_edges), strict=False)
+            nodes=io.StringIO(self.small_tree_ex_nodes),
+            edges=io.StringIO(self.small_tree_ex_edges), strict=False)
         tables = ts.dump_tables()
         nodes = tables.nodes
         flags = nodes.flags
@@ -2766,8 +2756,8 @@ def test_small_tree_internal_samples(self):
 
     def test_small_tree_linear_samples(self):
         ts = tskit.load_text(
-            nodes=six.StringIO(self.small_tree_ex_nodes),
-            edges=six.StringIO(self.small_tree_ex_edges), strict=False)
+            nodes=io.StringIO(self.small_tree_ex_nodes),
+            edges=io.StringIO(self.small_tree_ex_edges), strict=False)
         tables = ts.dump_tables()
         nodes = tables.nodes
         flags = nodes.flags
@@ -2788,8 +2778,8 @@ def test_small_tree_linear_samples(self):
 
     def test_small_tree_internal_and_external_samples(self):
         ts = tskit.load_text(
-            nodes=six.StringIO(self.small_tree_ex_nodes),
-            edges=six.StringIO(self.small_tree_ex_edges), strict=False)
+            nodes=io.StringIO(self.small_tree_ex_nodes),
+            edges=io.StringIO(self.small_tree_ex_edges), strict=False)
         tables = ts.dump_tables()
         nodes = tables.nodes
         flags = nodes.flags
@@ -2812,8 +2802,8 @@ def test_small_tree_internal_and_external_samples(self):
 
     def test_small_tree_mutations(self):
         ts = tskit.load_text(
-            nodes=six.StringIO(self.small_tree_ex_nodes),
-            edges=six.StringIO(self.small_tree_ex_edges), strict=False)
+            nodes=io.StringIO(self.small_tree_ex_nodes),
+            edges=io.StringIO(self.small_tree_ex_edges), strict=False)
         tables = ts.dump_tables()
         # Add some simple mutations here above the nodes we're keeping.
         tables.sites.add_row(position=0.25, ancestral_state="0")
@@ -2834,8 +2824,8 @@ def test_small_tree_mutations(self):
 
     def test_small_tree_filter_zero_mutations(self):
         ts = tskit.load_text(
-            nodes=six.StringIO(self.small_tree_ex_nodes),
-            edges=six.StringIO(self.small_tree_ex_edges), strict=False)
+            nodes=io.StringIO(self.small_tree_ex_nodes),
+            edges=io.StringIO(self.small_tree_ex_edges), strict=False)
         ts = tsutil.insert_branch_sites(ts)
         self.assertEqual(ts.num_sites, 8)
         self.assertEqual(ts.num_mutations, 8)
@@ -2848,8 +2838,8 @@ def test_small_tree_filter_zero_mutations(self):
 
     def test_small_tree_fixed_sites(self):
         ts = tskit.load_text(
-            nodes=six.StringIO(self.small_tree_ex_nodes),
-            edges=six.StringIO(self.small_tree_ex_edges), strict=False)
+            nodes=io.StringIO(self.small_tree_ex_nodes),
+            edges=io.StringIO(self.small_tree_ex_edges), strict=False)
         tables = ts.dump_tables()
         # Add some simple mutations that will be fixed after simplify
         tables.sites.add_row(position=0.25, ancestral_state="0")
@@ -2868,8 +2858,8 @@ def test_small_tree_fixed_sites(self):
 
     def test_small_tree_mutations_over_root(self):
         ts = tskit.load_text(
-            nodes=six.StringIO(self.small_tree_ex_nodes),
-            edges=six.StringIO(self.small_tree_ex_edges), strict=False)
+            nodes=io.StringIO(self.small_tree_ex_nodes),
+            edges=io.StringIO(self.small_tree_ex_edges), strict=False)
         tables = ts.dump_tables()
         tables.sites.add_row(position=0.25, ancestral_state="0")
         tables.mutations.add_row(site=0, node=8, derived_state="1")
@@ -2883,8 +2873,8 @@ def test_small_tree_mutations_over_root(self):
 
     def test_small_tree_recurrent_mutations(self):
         ts = tskit.load_text(
-            nodes=six.StringIO(self.small_tree_ex_nodes),
-            edges=six.StringIO(self.small_tree_ex_edges), strict=False)
+            nodes=io.StringIO(self.small_tree_ex_nodes),
+            edges=io.StringIO(self.small_tree_ex_edges), strict=False)
         tables = ts.dump_tables()
         # Add recurrent mutation on the root branches
         tables.sites.add_row(position=0.25, ancestral_state="0")
@@ -2901,8 +2891,8 @@ def test_small_tree_recurrent_mutations(self):
 
     def best_small_tree_back_mutations(self):
         ts = tskit.load_text(
-            nodes=six.StringIO(self.small_tree_ex_nodes),
-            edges=six.StringIO(self.small_tree_ex_edges), strict=False)
+            nodes=io.StringIO(self.small_tree_ex_nodes),
+            edges=io.StringIO(self.small_tree_ex_edges), strict=False)
         tables = ts.dump_tables()
         # Add a chain of mutations
         tables.sites.add_row(position=0.25, ancestral_state="0")
@@ -2935,13 +2925,13 @@ def best_small_tree_back_mutations(self):
         self.assertEqual(list(tss.haplotypes()), ["1", "0", "1"])
 
     def test_overlapping_unary_edges(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
         2       0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       2       2       0
         1       3       2       1
@@ -2957,13 +2947,13 @@ def test_overlapping_unary_edges(self):
             self.assertEqual(t.parent_dict, trees[t.index])
 
     def test_overlapping_unary_edges_internal_samples(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
         2       1           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       2       2       0
         1       3       2       1
@@ -2978,13 +2968,13 @@ def test_overlapping_unary_edges_internal_samples(self):
         self.assertEqual(list(node_map), [0, 1, 2])
 
     def test_isolated_samples(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           1
         2       1           2
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         """)
         ts = tskit.load_text(nodes, edges, sequence_length=1, strict=False)
@@ -2997,7 +2987,7 @@ def test_isolated_samples(self):
         self.assertEqual(list(node_map), [0, 1, 2])
 
     def test_internal_samples(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   population      time
         0       1       -1              1.00000000000000
         1       0       -1              1.00000000000000
@@ -3007,7 +2997,7 @@ def test_internal_samples(self):
         5       1       -1              0.00000000000000
 
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         id      left            right           parent  child
         0       0.62185118      1.00000000      1       5
         1       0.00000000      0.62185118      2       5
@@ -3030,20 +3020,20 @@ def test_internal_samples(self):
             self.assertEqual(t.parent_dict, trees[t.index])
 
     def test_many_mutations_over_single_sample_ancestral_state(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       1       0
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         position    ancestral_state
         0           0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state   parent
         0       0       1               -1
         0       0       0               0
@@ -3060,20 +3050,20 @@ def test_many_mutations_over_single_sample_ancestral_state(self):
         self.assertEqual(list(tss.haplotypes()), ["0"])
 
     def test_many_mutations_over_single_sample_derived_state(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       0           1
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0       1       1       0
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         position    ancestral_state
         0           0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state   parent
         0       0       1               -1
         0       0       0               0
@@ -3311,7 +3301,7 @@ def verify_parents(self, ts):
         self.assertTrue(np.array_equal(parent, tables.mutations.parent))
 
     def test_example(self):
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       0           2.0
         1       0           1.0
@@ -3319,7 +3309,7 @@ def test_example(self):
         3       1           0
         4       1           0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.0    0.5   2  3
         0.0    0.8   2  4
@@ -3328,13 +3318,13 @@ def test_example(self):
         0.0    1.0   0  2
         0.8    1.0   0  4
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         position    ancestral_state
         0.1     0
         0.5     0
         0.9     0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state   parent
         0       1       1               -1
         0       2       1               -1
diff --git a/python/tests/test_tree_stats.py b/python/tests/test_tree_stats.py
index f37269304f..7c94ef0c6c 100644
--- a/python/tests/test_tree_stats.py
+++ b/python/tests/test_tree_stats.py
@@ -23,17 +23,13 @@
 """
 Test cases for generalized statistic computation.
 """
-from __future__ import print_function
-from __future__ import division
-
-
+import io
 import unittest
 import random
 
 import numpy as np
 import numpy.testing as nt
 
-import six
 import msprime
 
 import tskit
@@ -144,7 +140,7 @@ def Y2(self, X, Y, begin=0.0, end=None):
             this_length = min(end, tr.interval[1]) - max(begin, tr.interval[0])
             for x in X:
                 for y in Y:
-                    for z in set(Y) - set([y]):
+                    for z in set(Y) - {y}:
                         xy_mrca = tr.mrca(x, y)
                         xz_mrca = tr.mrca(x, z)
                         yz_mrca = tr.mrca(y, z)
@@ -176,8 +172,8 @@ def Y1(self, X, begin=0.0, end=None):
                 break
             this_length = min(end, tr.interval[1]) - max(begin, tr.interval[0])
             for x in X:
-                for y in set(X) - set([x]):
-                    for z in set(X) - set([x, y]):
+                for y in set(X) - {x}:
+                    for z in set(X) - {x, y}:
                         xy_mrca = tr.mrca(x, y)
                         xz_mrca = tr.mrca(x, z)
                         yz_mrca = tr.mrca(y, z)
@@ -239,7 +235,7 @@ def f3(self, A, B, C, begin=0.0, end=None):
             SS = 0
             for a in A:
                 for b in B:
-                    for c in set(A) - set([a]):
+                    for c in set(A) - {a}:
                         for d in C:
                             SS += path_length(tr, tr.mrca(a, c), tr.mrca(b, d))
                             SS -= path_length(tr, tr.mrca(a, d), tr.mrca(b, c))
@@ -264,8 +260,8 @@ def f2(self, A, B, begin=0.0, end=None):
             SS = 0
             for a in A:
                 for b in B:
-                    for c in set(A) - set([a]):
-                        for d in set(B) - set([b]):
+                    for c in set(A) - {a}:
+                        for d in set(B) - {b}:
                             SS += path_length(tr, tr.mrca(a, c), tr.mrca(b, d))
                             SS -= path_length(tr, tr.mrca(a, d), tr.mrca(b, c))
             S += SS * this_length
@@ -400,7 +396,7 @@ def Y2(self, X, Y, begin=0.0, end=None):
             if (site_positions[k] >= begin) and (site_positions[k] < end):
                 for x in X:
                     for y in Y:
-                        for z in set(Y) - set([y]):
+                        for z in set(Y) - {y}:
                             if ((haps[x][k] != haps[y][k])
                                and (haps[x][k] != haps[z][k])):
                                 # x|yz
@@ -416,8 +412,8 @@ def Y1(self, X, begin=0.0, end=None):
         for k in range(self.tree_sequence.num_sites):
             if (site_positions[k] >= begin) and (site_positions[k] < end):
                 for x in X:
-                    for y in set(X) - set([x]):
-                        for z in set(X) - set([x, y]):
+                    for y in set(X) - {x}:
+                        for z in set(X) - {x, y}:
                             if ((haps[x][k] != haps[y][k])
                                and (haps[x][k] != haps[z][k])):
                                 # x|yz
@@ -464,7 +460,7 @@ def f3(self, A, B, C, begin=0.0, end=None):
             if (site_positions[k] >= begin) and (site_positions[k] < end):
                 for a in A:
                     for b in B:
-                        for c in set(A) - set([a]):
+                        for c in set(A) - {a}:
                             for d in C:
                                 if ((haps[a][k] == haps[c][k])
                                    and (haps[a][k] != haps[d][k])
@@ -491,8 +487,8 @@ def f2(self, A, B, begin=0.0, end=None):
             if (site_positions[k] >= begin) and (site_positions[k] < end):
                 for a in A:
                     for b in B:
-                        for c in set(A) - set([a]):
-                            for d in set(B) - set([b]):
+                        for c in set(A) - {a}:
+                            for d in set(B) - {b}:
                                 if ((haps[a][k] == haps[c][k])
                                    and (haps[a][k] != haps[d][k])
                                    and (haps[a][k] != haps[b][k])):
@@ -836,7 +832,7 @@ def test_case_1(self):
         branch_true_Y = 0.2*(1 + 0.5) + 0.6*(0.4) + 0.2*(0.7+0.2)
         site_true_Y = 3 + 0 + 1
 
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -846,7 +842,7 @@ def test_case_1(self):
         5       0           0.7
         6       0           1.0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.2     0.8     3       0,2
         0.0     0.2     4       1,2
@@ -855,7 +851,7 @@ def test_case_1(self):
         0.8     1.0     5       0,4
         0.0     0.2     6       0,4
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         id  position    ancestral_state
         0   0.05        0
         1   0.1         0
@@ -868,7 +864,7 @@ def test_case_1(self):
         8   0.95        0
         9   0.951       0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state
         0       4       1
         1       0       1
@@ -947,23 +943,23 @@ def f(x):
     def test_case_odds_and_ends(self):
         # Tests having (a) the first site after the first window, and
         # (b) no samples having the ancestral state.
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
         2       0           0.5
         3       0           1.0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.0     0.5     2       0,1
         0.5     1.0     3       0,1
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         id  position    ancestral_state
         0   0.65        0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state   parent
         0       0       1               -1
         0       1       2               -1
@@ -997,7 +993,7 @@ def test_case_recurrent_muts(self):
         #       0     2       0        1   0   1       0   2       3
         site_true_Y = 0 + 1 + 1
 
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         id      is_sample   time
         0       1           0
         1       1           0
@@ -1007,7 +1003,7 @@ def test_case_recurrent_muts(self):
         5       0           0.7
         6       0           1.0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.2     0.8     3       0,2
         0.0     0.2     4       1,2
@@ -1016,13 +1012,13 @@ def test_case_recurrent_muts(self):
         0.8     1.0     5       0,4
         0.0     0.2     6       0,4
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         id  position    ancestral_state
         0   0.05        0
         1   0.3         0
         2   0.9         0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state   parent
         0       0       1               -1
         0       0       2               0
@@ -1098,7 +1094,7 @@ def test_case_2(self):
         # Y(0;1, 2)
         site_true_Y = 1
 
-        nodes = six.StringIO("""\
+        nodes = io.StringIO("""\
         is_sample       time    population
         1       0.000000        0
         1       0.000000        0
@@ -1112,7 +1108,7 @@ def test_case_2(self):
         0       1.000000        0
         0       1.000000        0
         """)
-        edges = six.StringIO("""\
+        edges = io.StringIO("""\
         left    right   parent  child
         0.500000        1.000000        10      1
         0.000000        0.400000        10      2
@@ -1133,14 +1129,14 @@ def test_case_2(self):
         0.100000        0.900000        3       4,5
         0.000000        0.100000        3       4,5,7
         """)
-        sites = six.StringIO("""\
+        sites = io.StringIO("""\
         id  position    ancestral_state
         0   0.0         0
         1   0.55        0
         2   0.75        0
         3   0.85        0
         """)
-        mutations = six.StringIO("""\
+        mutations = io.StringIO("""\
         site    node    derived_state   parent
         0       0       1               -1
         0       10      1               -1
diff --git a/python/tests/test_vcf.py b/python/tests/test_vcf.py
index e391651b61..b62bbb7e94 100644
--- a/python/tests/test_vcf.py
+++ b/python/tests/test_vcf.py
@@ -23,9 +23,6 @@
 """
 Test cases for VCF output in tskit.
 """
-from __future__ import print_function
-from __future__ import division
-
 import collections
 import math
 import os
diff --git a/python/tests/test_wright_fisher.py b/python/tests/test_wright_fisher.py
index 8f079695b2..d2f5698377 100644
--- a/python/tests/test_wright_fisher.py
+++ b/python/tests/test_wright_fisher.py
@@ -23,9 +23,6 @@
 """
 Test various functions using messy tables output by a forwards-time simulator.
 """
-from __future__ import print_function
-from __future__ import division
-
 import itertools
 import random
 import unittest
diff --git a/python/tests/tsutil.py b/python/tests/tsutil.py
index 33d1ed1cd5..2c2faf15c2 100644
--- a/python/tests/tsutil.py
+++ b/python/tests/tsutil.py
@@ -23,10 +23,6 @@
 """
 A collection of utilities to edit and construct tree sequences.
 """
-from __future__ import print_function
-from __future__ import unicode_literals
-from __future__ import division
-
 import json
 import random
 
diff --git a/python/tskit/__init__.py b/python/tskit/__init__.py
index 4a188ab54f..41d94b2ab9 100644
--- a/python/tskit/__init__.py
+++ b/python/tskit/__init__.py
@@ -20,9 +20,6 @@
 # OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
 # SOFTWARE.
 
-from __future__ import print_function
-from __future__ import division
-
 import _tskit
 FORWARD = _tskit.FORWARD
 REVERSE = _tskit.REVERSE
diff --git a/python/tskit/__main__.py b/python/tskit/__main__.py
index a0f7cb8179..b07a50eba2 100644
--- a/python/tskit/__main__.py
+++ b/python/tskit/__main__.py
@@ -1,7 +1,3 @@
-from __future__ import print_function
-from __future__ import division
-
-
 def main():
     print("Not implemented yet; please be patient!!")
 
diff --git a/python/tskit/drawing.py b/python/tskit/drawing.py
index 43343c029d..1b7b98d427 100644
--- a/python/tskit/drawing.py
+++ b/python/tskit/drawing.py
@@ -23,12 +23,8 @@
 """
 Module responsible for visualisations.
 """
-from __future__ import division
-from __future__ import print_function
-
 import array
 import collections
-import sys
 
 try:
     import svgwrite
@@ -36,8 +32,6 @@
 except ImportError:  # pragma: no cover
     _svgwrite_imported = False
 
-IS_PY2 = sys.version_info[0] < 3
-
 NULL_NODE = -1
 
 
@@ -64,8 +58,6 @@ def draw_tree(
     elif fmt == "ascii":
         cls = AsciiTreeDrawer
     elif fmt == "unicode":
-        if IS_PY2:
-            raise ValueError("Unicode tree drawing not supported on Python 2")
         cls = UnicodeTreeDrawer
 
     # We can't draw trees with zero roots.
@@ -294,7 +286,7 @@ def _convert_text(self, text):
 
     def _assign_coordinates(self):
         # Get the age of each node and rank them.
-        times = set(self._tree.time(u) for u in self._tree.nodes())
+        times = {self._tree.time(u) for u in self._tree.nodes()}
         depth = {t: 2 * j for j, t in enumerate(sorted(times, reverse=True))}
         for u in self._tree.nodes():
             self._y_coords[u] = depth[self._tree.time(u)]
@@ -387,9 +379,7 @@ def _convert_text(self, text):
         return array.array(self.array_type, text.encode())
 
     def draw(self):
-        s = self._draw().tostring()
-        if not IS_PY2:
-            s = s.decode()
+        s = self._draw().tostring().decode()
         return s
 
 
diff --git a/python/tskit/exceptions.py b/python/tskit/exceptions.py
index a4b74dfaac..80b957ffb6 100644
--- a/python/tskit/exceptions.py
+++ b/python/tskit/exceptions.py
@@ -36,19 +36,15 @@
 
 # TODO finalise this when working out the docs structure for tskit on rtd.
 
-try:
-    TskitException.__doc__ = "Superclass of all exceptions defined in tskit."
-    LibraryError.__doc__ = "Generic low-level error raised by the C library."
-    FileFormatError.__doc__ = "An error was detected in the file format."
-    VersionTooNewError.__doc__ = """
-    The version of the file is too new and cannot be read by the library.
-    """
-    VersionTooOldError.__doc__ = """
-    The version of the file is too old and cannot be read by the library.
-    """
-except AttributeError:  # pragma: no cover
-    # Python2 throws attribute error. Ignore.
-    pass
+TskitException.__doc__ = "Superclass of all exceptions defined in tskit."
+LibraryError.__doc__ = "Generic low-level error raised by the C library."
+FileFormatError.__doc__ = "An error was detected in the file format."
+VersionTooNewError.__doc__ = """
+The version of the file is too new and cannot be read by the library.
+"""
+VersionTooOldError.__doc__ = """
+The version of the file is too old and cannot be read by the library.
+"""
 
 
 class DuplicatePositionsError(TskitException):
diff --git a/python/tskit/formats.py b/python/tskit/formats.py
index 6d8403fae5..81f0e42f54 100644
--- a/python/tskit/formats.py
+++ b/python/tskit/formats.py
@@ -24,9 +24,6 @@
 Module responsible for converting tree sequence files from older
 formats.
 """
-from __future__ import division
-from __future__ import print_function
-
 import datetime
 import json
 import logging
diff --git a/python/tskit/provenance.py b/python/tskit/provenance.py
index a99ebb9f3a..9dbab6bc65 100644
--- a/python/tskit/provenance.py
+++ b/python/tskit/provenance.py
@@ -24,9 +24,6 @@
 Common provenance methods used to determine the state and versions
 of various dependencies and the OS.
 """
-from __future__ import print_function
-from __future__ import division
-
 import platform
 import json
 import os.path
diff --git a/python/tskit/stats.py b/python/tskit/stats.py
index 63acb153c8..5d415b323d 100644
--- a/python/tskit/stats.py
+++ b/python/tskit/stats.py
@@ -22,9 +22,6 @@
 """
 Module responsible for computing various statistics on tree sequences.
 """
-from __future__ import division
-from __future__ import print_function
-
 import threading
 import struct
 import sys
diff --git a/python/tskit/tables.py b/python/tskit/tables.py
index 9d82660b78..e95ffaed8a 100644
--- a/python/tskit/tables.py
+++ b/python/tskit/tables.py
@@ -23,16 +23,12 @@
 """
 Tree sequence IO via the tables API.
 """
-from __future__ import division
-from __future__ import print_function
-
 import base64
 import collections
 import datetime
 import warnings
 
 import numpy as np
-from six.moves import copyreg
 
 import _tskit
 # This circular import is ugly but it seems hard to avoid it since table collection
@@ -149,8 +145,13 @@ def truncate(self, num_rows):
         """
         return self.ll_table.truncate(num_rows)
 
+    # Picle support
+    def __getstate__(self):
+        return self.asdict()
+
     # Unpickle support
     def __setstate__(self, state):
+        self.__init__()
         self.set_columns(**state)
 
     def asdict(self):
@@ -193,7 +194,7 @@ class IndividualTable(BaseTable):
     def __init__(self, max_rows_increment=0, ll_table=None):
         if ll_table is None:
             ll_table = _tskit.IndividualTable(max_rows_increment=max_rows_increment)
-        super(IndividualTable, self).__init__(ll_table, IndividualTableRow)
+        super().__init__(ll_table, IndividualTableRow)
 
     @property
     def flags(self):
@@ -336,11 +337,6 @@ def asdict(self):
         }
 
 
-# Pickle support. See copyreg registration for this function below.
-def _pickle_individual_table(table):
-    return IndividualTable, tuple(), table.asdict()
-
-
 class NodeTable(BaseTable):
     """
     A table defining the nodes in a tree sequence. See the
@@ -373,7 +369,7 @@ class NodeTable(BaseTable):
     def __init__(self, max_rows_increment=0, ll_table=None):
         if ll_table is None:
             ll_table = _tskit.NodeTable(max_rows_increment=max_rows_increment)
-        super(NodeTable, self).__init__(ll_table, NodeTableRow)
+        super().__init__(ll_table, NodeTableRow)
 
     @property
     def time(self):
@@ -532,11 +528,6 @@ def asdict(self):
         }
 
 
-# Pickle support. See copyreg registration for this function below.
-def _pickle_node_table(table):
-    return NodeTable, tuple(), table.asdict()
-
-
 class EdgeTable(BaseTable):
     """
     A table defining the edges in a tree sequence. See the
@@ -563,7 +554,7 @@ class EdgeTable(BaseTable):
     def __init__(self, max_rows_increment=0, ll_table=None):
         if ll_table is None:
             ll_table = _tskit.EdgeTable(max_rows_increment=max_rows_increment)
-        super(EdgeTable, self).__init__(ll_table, EdgeTableRow)
+        super().__init__(ll_table, EdgeTableRow)
 
     @property
     def left(self):
@@ -666,11 +657,6 @@ def asdict(self):
         }
 
 
-# Pickle support. See copyreg registration for this function below.
-def _edge_table_pickle(table):
-    return EdgeTable, tuple(), table.asdict()
-
-
 class MigrationTable(BaseTable):
     """
     A table defining the migrations in a tree sequence. See the
@@ -702,7 +688,7 @@ class MigrationTable(BaseTable):
     def __init__(self, max_rows_increment=0, ll_table=None):
         if ll_table is None:
             ll_table = _tskit.MigrationTable(max_rows_increment=max_rows_increment)
-        super(MigrationTable, self).__init__(ll_table, MigrationTableRow)
+        super().__init__(ll_table, MigrationTableRow)
 
     @property
     def left(self):
@@ -830,11 +816,6 @@ def asdict(self):
         }
 
 
-# Pickle support. See copyreg registration for this function below.
-def _migration_table_pickle(table):
-    return MigrationTable, tuple(), table.asdict()
-
-
 class SiteTable(BaseTable):
     """
     A table defining the sites in a tree sequence. See the
@@ -868,7 +849,7 @@ class SiteTable(BaseTable):
     def __init__(self, max_rows_increment=0, ll_table=None):
         if ll_table is None:
             ll_table = _tskit.SiteTable(max_rows_increment=max_rows_increment)
-        super(SiteTable, self).__init__(ll_table, SiteTableRow)
+        super().__init__(ll_table, SiteTableRow)
 
     @property
     def position(self):
@@ -1017,11 +998,6 @@ def asdict(self):
         }
 
 
-# Pickle support. See copyreg registration for this function below.
-def _site_table_pickle(table):
-    return SiteTable, tuple(), table.asdict()
-
-
 class MutationTable(BaseTable):
     """
     A table defining the mutations in a tree sequence. See the
@@ -1059,7 +1035,7 @@ class MutationTable(BaseTable):
     def __init__(self, max_rows_increment=0, ll_table=None):
         if ll_table is None:
             ll_table = _tskit.MutationTable(max_rows_increment=max_rows_increment)
-        super(MutationTable, self).__init__(ll_table, MutationTableRow)
+        super().__init__(ll_table, MutationTableRow)
 
     @property
     def site(self):
@@ -1232,11 +1208,6 @@ def asdict(self):
         }
 
 
-# Pickle support. See copyreg registration for this function below.
-def _mutation_table_pickle(table):
-    return MutationTable, tuple(), table.asdict()
-
-
 class PopulationTable(BaseTable):
     """
     A table defining the populations referred to in a tree sequence.
@@ -1262,7 +1233,7 @@ class PopulationTable(BaseTable):
     def __init__(self, max_rows_increment=0, ll_table=None):
         if ll_table is None:
             ll_table = _tskit.PopulationTable(max_rows_increment=max_rows_increment)
-        super(PopulationTable, self).__init__(ll_table, PopulationTableRow)
+        super().__init__(ll_table, PopulationTableRow)
 
     @property
     def metadata(self):
@@ -1317,11 +1288,6 @@ def asdict(self):
         }
 
 
-# Pickle support. See copyreg registration for this function below.
-def _population_table_pickle(table):
-    return PopulationTable, tuple(), table.asdict()
-
-
 class ProvenanceTable(BaseTable):
     """
     A table recording the provenance (i.e., history) of this table, so that the
@@ -1350,7 +1316,7 @@ class ProvenanceTable(BaseTable):
     def __init__(self, max_rows_increment=0, ll_table=None):
         if ll_table is None:
             ll_table = _tskit.ProvenanceTable(max_rows_increment=max_rows_increment)
-        super(ProvenanceTable, self).__init__(ll_table, ProvenanceTableRow)
+        super().__init__(ll_table, ProvenanceTableRow)
 
     @property
     def record(self):
@@ -1408,14 +1374,6 @@ def __str__(self):
             ret += "{}\t{}\t{}\n".format(j, timestamp[j], record[j])
         return ret[:-1]
 
-    # Unpickle support
-    def __setstate__(self, state):
-        self.set_columns(
-            timestamp=state["timestamp"],
-            timestamp_offset=state["timestamp_offset"],
-            record=state["record"],
-            record_offset=state["record_offset"])
-
     def copy(self):
         """
         Returns a deep copy of this table.
@@ -1437,11 +1395,6 @@ def asdict(self):
         }
 
 
-# Pickle support. See copyreg registration for this function below.
-def _provenance_table_pickle(table):
-    return ProvenanceTable, tuple(), table.asdict()
-
-
 class TableCollection(object):
     """
     A collection of mutable tables defining a tree sequence. See the
@@ -1580,6 +1533,9 @@ def __eq__(self, other):
     def __ne__(self, other):
         return not self.__eq__(other)
 
+    def __getstate__(self):
+        return self.asdict()
+
     # Unpickle support
     def __setstate__(self, state):
         self.__init__(state["sequence_length"])
@@ -1776,26 +1732,6 @@ def drop_index(self):
         self.ll_tables.drop_index()
 
 
-# Pickle support. See copyreg registration for this function below.
-def _table_collection_pickle(tables):
-    return TableCollection, tuple(), tables.asdict()
-
-
-# Pickle support for the various tables. We are forced to use copyreg.pickle
-# here to support Python 2. For Python 3, we can just use the __setstate__.
-# It would be cleaner to attach the pickle_*_table functions to the classes
-# themselves, but this causes issues with Mocking on readthedocs. Sigh.
-copyreg.pickle(IndividualTable, _pickle_individual_table)
-copyreg.pickle(NodeTable, _pickle_node_table)
-copyreg.pickle(EdgeTable, _edge_table_pickle)
-copyreg.pickle(MigrationTable, _migration_table_pickle)
-copyreg.pickle(SiteTable, _site_table_pickle)
-copyreg.pickle(MutationTable, _mutation_table_pickle)
-copyreg.pickle(PopulationTable, _population_table_pickle)
-copyreg.pickle(ProvenanceTable, _provenance_table_pickle)
-copyreg.pickle(TableCollection, _table_collection_pickle)
-
-
 #############################################
 # Table functions.
 #############################################
diff --git a/python/tskit/trees.py b/python/tskit/trees.py
index f5aa0e00a5..43650b24f6 100644
--- a/python/tskit/trees.py
+++ b/python/tskit/trees.py
@@ -25,21 +25,13 @@
 """
 Module responsible for managing trees and tree sequences.
 """
-from __future__ import division
-from __future__ import print_function
-
 import collections
 import itertools
 import json
-import sys
 import base64
 import warnings
 import functools
-try:
-    import concurrent.futures
-except ImportError:  # pragma: no cover
-    # We're on Python2; any attempts to use futures are dealt with below.
-    pass
+import concurrent.futures
 
 import numpy as np
 
@@ -54,9 +46,6 @@
 from _tskit import NULL
 
 
-IS_PY2 = sys.version_info[0] < 3
-
-
 CoalescenceRecord = collections.namedtuple(
     "CoalescenceRecord",
     ["left", "right", "node", "children", "time", "population"])
@@ -1306,7 +1295,7 @@ def __build_newick(self, node, precision, node_labels):
         """
         label = node_labels.get(node, "")
         if self.is_leaf(node):
-            s = "{0}".format(label)
+            s = "{}".format(label)
         else:
             s = "("
             for child in self.children(node):
@@ -1350,8 +1339,7 @@ def newick(self, precision=14, root=None, node_labels=None):
             root = self.root
         if node_labels is None:
             s = self._ll_tree.get_newick(precision=precision, root=root)
-            if not IS_PY2:
-                s = s.decode()
+            s = s.decode()
         else:
             return self.__build_newick(root, precision, node_labels) + ";"
         return s
@@ -1833,10 +1821,6 @@ def __reversed__(self):
         self.forward = False
         return self
 
-    def next(self):
-        # This is for Python2 support only.
-        return self.__next__()
-
     def __next__(self):
         if self.forward:
             self.more_trees = self.more_trees and self.tree.next()
@@ -2541,11 +2525,8 @@ def haplotypes(self):
             multiletter alleles.
         """
         hapgen = _tskit.HaplotypeGenerator(self._ll_tree_sequence)
-        j = 0
-        # Would use range here except for Python 2.
-        while j < self.num_samples:
+        for j in range(self.num_samples):
             yield hapgen.get_haplotype(j)
-            j += 1
 
     # Samples is experimental for now, so we don't document it.
     def variants(self, as_bytes=False, samples=None):
@@ -2664,8 +2645,6 @@ def genealogical_nearest_neighbours(self, focal, reference_sets, num_threads=0):
             return self._ll_tree_sequence.genealogical_nearest_neighbours(
                 focal, reference_sets)
         else:
-            if IS_PY2:  # pragma: no cover
-                raise ValueError("Threads not supported on Python 2.")
             worker = functools.partial(
                 self._ll_tree_sequence.genealogical_nearest_neighbours,
                 reference_sets=reference_sets)