diff --git a/.gitignore b/.gitignore
index 5cbbeeb..565fdcb 100755
--- a/.gitignore
+++ b/.gitignore
@@ -39,3 +39,4 @@ output/*/index.html
 
 # Sphinx
 docs/_build
+tutorials/1_basic_popgen_setup/*all_controls_ipu
diff --git a/.travis.yml b/.travis.yml
index 5bb737e..60b36f9 100755
--- a/.travis.yml
+++ b/.travis.yml
@@ -4,33 +4,38 @@ language: python
 
 python:
   - "2.7"
+  - "3.5"
+  - "3.6"
 
-#before_install:
-#  - sudo apt-get install liblapack-dev gfortran
-
-# command to install dependencies, e.g. pip install -r requirements.txt --use-mirrors
-# install:
-#  - travis_wait pip install -r requirements.txt
-
-
-before_install:
-  - wget http://repo.continuum.io/miniconda/Miniconda-latest-Linux-x86_64.sh -O miniconda.sh
-  - chmod +x miniconda.sh
-  - ./miniconda.sh -b
-  - export PATH=/home/travis/miniconda/bin:$PATH
-  # Update conda itself
-  - conda update --yes conda
-  # The next couple lines fix a crash with multiprocessing
-  # on Travis and are not specific to using Miniconda
-  - sudo rm -rf /dev/shm
-  - sudo ln -s /run/shm /dev/shm
 
 install:
-  - conda create --yes -n popgen_venv python=$TRAVIS_PYTHON_VERSION pip numpy scipy pandas PyYAML
+  - sudo apt-get update
+  # We do this conditionally because it saves us some downloading if the
+  # version is the same.
+  - if [[ "$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
+  - bash miniconda.sh -b -p $HOME/miniconda
+  - export PATH="$HOME/miniconda/bin:$PATH"
+  - hash -r
+  - conda config --set always_yes yes --set changeps1 no
+  - conda update -q conda
+  # Useful for debugging any issues with conda
+  - conda info -a
+
+  - conda create -q -n popgen_venv python=$TRAVIS_PYTHON_VERSION pip numpy scipy pandas PyYAML pytest
   - source activate popgen_venv
-  - pip install pytest
   - python setup.py install
 
 # command to run tests, e.g. python setup.py test
 script:
   - py.test
+
+
+
+
+
+
+
diff --git a/AUTHORS.rst b/AUTHORS.rst
index 6478b32..37def8a 100755
--- a/AUTHORS.rst
+++ b/AUTHORS.rst
@@ -11,3 +11,5 @@ Contributors
 ------------
 
 * Daehyun You
+* Jeffrey Newman <jnewman@camsys.com>
+
diff --git a/HISTORY.rst b/HISTORY.rst
index 4aa06c3..3015d3a 100755
--- a/HISTORY.rst
+++ b/HISTORY.rst
@@ -3,7 +3,13 @@
 History
 -------
 
+2.1.b1 (March 23, 2017)
++++++++++++++++++++++++
+
+* Updated for Python 3.5+
+
+
 2.0.b1 (June 4, 2015)
-++++++++++++++++++
++++++++++++++++++++++
 
 * This is the first non-production release of PopGen 2.0
diff --git a/README.rst b/README.rst
index 380698e..ec4fe7e 100755
--- a/README.rst
+++ b/README.rst
@@ -1,6 +1,6 @@
-=============================
+==================================
 Synthetic Population Generator 2.0
-=============================
+==================================
 
 .. image:: https://badge.fury.io/py/popgen.png
     :target: http://badge.fury.io/py/popgen
@@ -12,15 +12,6 @@ Synthetic Population Generator 2.0
     :target: https://pypi.python.org/pypi/popgen
 
 
-Synthetic Population Generator 2.0
-
-
-Features
---------
 
-For PopGen 2.1
-* Automatic data validator
-* Census data importer
+This fork contains code that has been modified to be compatible with Python 2.7 and 3.5 and up.
 
-For PopGen 2.2
-* Output Data Analyzer with charting and table generation capabilities
diff --git a/bin/popgen_run b/bin/popgen_run
index bd1587d..018622c 100644
--- a/bin/popgen_run
+++ b/bin/popgen_run
@@ -18,7 +18,7 @@ def run():
     """
     args = sys.argv[1:]
 
-    print args
+    print (args)
 
     if len(args) != 1:
         raise ArgumentsError(
@@ -35,6 +35,6 @@ def run():
     p_obj = Project(config_loc=args[0])
     p_obj.load_project()
     p_obj.run_scenarios()
-    print "Time it took: %.4f" % (time.time() - t)
+    print ("Time it took: %.4f" % (time.time() - t))
 
 run()
diff --git a/popgen/__init__.py b/popgen/__init__.py
index 5ce17c6..f8fd35e 100755
--- a/popgen/__init__.py
+++ b/popgen/__init__.py
@@ -1 +1,3 @@
 from .project import Project
+
+__version__ = '2.1.b2'  # now for Python 3!
diff --git a/popgen/config.py b/popgen/config.py
old mode 100644
new mode 100755
index bfeeae4..b7b3bee
--- a/popgen/config.py
+++ b/popgen/config.py
@@ -1,3 +1,4 @@
+from __future__ import print_function
 import yaml
 
 
@@ -16,7 +17,7 @@ def wrap_config_value(value):
     then attribute access x.attribute1.attribute2 is used to access "Value".
     Also, x.attribute can be used to access the dictionary {attribute: 'value'}
     """
-    if isinstance(value, basestring):
+    if isinstance(value, str):
         return value
 
     try:
@@ -45,7 +46,7 @@ def __getitem__(self, key):
     def return_value(self, key):
         try:
             value = self._data[key]
-        except KeyError, e:
+        except KeyError as e:
             raise ConfigError(
                 "Key - %s doesn't exist in the YAML configuration" % key)
         return value
@@ -78,5 +79,5 @@ def write_to_file(self, filepath):
 
     config_dict = yaml.load(yaml_f)
     config_obj = Config(config_dict)
-    print config_obj.project.name
-    print config_obj["project"]["name"]
+    print (config_obj.project.name)
+    print (config_obj["project"]["name"])
diff --git a/popgen/data.py b/popgen/data.py
index a532c98..2f20f22 100755
--- a/popgen/data.py
+++ b/popgen/data.py
@@ -1,9 +1,10 @@
+from __future__ import print_function
 import os
 
 import pandas as pd
 import numpy as np
 
-from config import ConfigError
+from .config import ConfigError
 
 
 class DB(object):
@@ -89,8 +90,8 @@ def enumerate_geo_ids_for_scenario(self, scenario_config):
             self.geo_ids = []
             for region_id in self.region_ids:
                 self.geo_ids += self.get_geo_ids_for_region(region_id)
-        except ConfigError, e:
-            print "KeyError", e
+        except ConfigError as e:
+            print ("KeyError", e)
             self.geo_ids = self.geo_ids_all
             # self.sample_geo_ids = self.sample_geo_ids_all
             self.region_ids = self.region_ids_all
diff --git a/popgen/draw.py b/popgen/draw.py
index 0e394c2..7125dc3 100755
--- a/popgen/draw.py
+++ b/popgen/draw.py
@@ -1,3 +1,4 @@
+from __future__ import print_function
 import pandas as pd
 import numpy as np
 from scipy import stats
@@ -86,11 +87,11 @@ def _return_cumulative_probability(self, geo_sample_weights):
 
     def _pick_households(self, geo_id_frequencies, geo_cumulative_weights):
         last = 0
-        rand_numbers = np.random.random(geo_id_frequencies.sum())
+        rand_numbers = np.random.random(int(geo_id_frequencies.sum()))
         list_rows_syn_subpop = []
         for column in self.geo_frequencies.columns.values:
             rows = self.geo_row_idx[column]
-            column_frequency = geo_id_frequencies[column]
+            column_frequency = int(geo_id_frequencies[column])
             column_bins = np.searchsorted(geo_cumulative_weights[column],
                                           (rand_numbers[
                                            last:last + column_frequency]),
diff --git a/popgen/ipf.py b/popgen/ipf.py
index 5bce334..63b97d4 100755
--- a/popgen/ipf.py
+++ b/popgen/ipf.py
@@ -1,3 +1,4 @@
+from __future__ import print_function
 import copy
 
 import numpy as np
@@ -62,7 +63,7 @@ def __init__(self, seed_all, seed, idx, marginals, ipf_config,
     def run_ipf(self):
         self.frequencies = self._correct_zero_cell_issue()
         # self.frequencies = self.seed["frequency"].values
-        for c_iter in xrange(self.ipf_iters):
+        for c_iter in range(self.ipf_iters):
             # print "Iter:", c_iter
             self._adjust_cell_frequencies()
             # Checks for convergence every 5 iterations
@@ -273,7 +274,7 @@ def _run_ipf_all_geos(self, entity, seed_geo, seed_all, row_idx, marginals,
             sample_geo_id = geo_corr_to_sample.loc[geo_id,
                                                    self.sample_geo_name]
             if isinstance(sample_geo_id, pd.Series):
-                seed_geo_levels_list = range(len(seed_geo.index.names))
+                seed_geo_levels_list = tuple(range(len(seed_geo.index.names)))
                 seed_for_geo_id = (seed_geo.loc[sample_geo_id.tolist()]
                                    .sum(level=seed_geo_levels_list[1:]))
                 # print (seed_geo.loc[sample_geo_id.tolist()])
@@ -340,7 +341,7 @@ def _get_stacked_constraints(self, constraints_list):
 
     def _get_columns_constraints_dict(self, constraints_dict):
         columns_constraints_dict = {}
-        for entity, constraints in constraints_dict.iteritems():
+        for entity, constraints in constraints_dict.items():
             columns_constraints_dict[entity] = (constraints
                                                 .index.values.tolist())
         # print columns_constraints_dict
diff --git a/popgen/output.py b/popgen/output.py
old mode 100644
new mode 100755
index cc796ab..428f891
--- a/popgen/output.py
+++ b/popgen/output.py
@@ -1,9 +1,12 @@
+from __future__ import print_function
 import os
 import time
+from warnings import warn
 
 import pandas as pd
 import numpy as np
 
+from .config import ConfigError
 
 class Syn_Population(object):
     def __init__(
@@ -125,16 +128,16 @@ def _get_stacked_sample(self, entities):
         sample_list = [self.db.sample[entity]
                        for entity in entities]
         stacked_sample = pd.concat(sample_list).fillna(0)
-        stacked_sample.sort(inplace=True)
+        stacked_sample.sort_index(inplace=True)
         return stacked_sample
 
     def add_records(self):
         geo_id_rows_syn_dict = self.draw_population_obj.geo_id_rows_syn_dict
-        for geo_id, geo_id_rows_syn in geo_id_rows_syn_dict.iteritems():
+        for geo_id, geo_id_rows_syn in geo_id_rows_syn_dict.items():
             geo_id_pop_syn = (
                 self._get_stacked_geo_for_geo_id(geo_id, geo_id_rows_syn))
             self.pop_rows_syn_dict[geo_id] = geo_id_pop_syn
-            self.pop_rows_syn_dict[geo_id][self.unique_id_in_geo_name] = (
+            self.pop_rows_syn_dict[geo_id][self.unique_id_in_geo_name] = tuple(
                 range(1, geo_id_rows_syn.shape[0]+1))
             # print self.pop_rows_syn_dict[geo_id].head()
             # raw_input()
@@ -166,14 +169,14 @@ def _get_housing_data_for_indexes(self, geo_id_pop_syn):
         housing_data = (
             geo_id_pop_syn.loc[:, [self.geo_name]].join(
                 self.housing_stacked_sample))
-        print "housing rows:", housing_data.shape
+        print ("housing rows:", housing_data.shape)
         return housing_data
 
     def _get_person_data_for_indexes(self, geo_id_pop_syn):
         person_data = (
             geo_id_pop_syn.loc[:, [self.geo_name]].join(
                 self.person_stacked_sample))
-        print "person rows:", person_data.shape
+        print ("person rows:", person_data.shape)
         return person_data
 
     def prepare_data(self):
@@ -191,8 +194,8 @@ def _stack_records(self):
         #    self.person_syn_dict.values(), copy=False)
         # print self.pop_syn_data["housing"].shape
         # print self.pop_syn_data["person"].shape
-        print "Time elapsed for stacking population is : %.4f" % (
-            time.time() - t)
+        print ("Time elapsed for stacking population is : %.4f" % (
+            time.time() - t))
 
     def _create_synthetic_population(self):
         t = time.time()
@@ -202,10 +205,10 @@ def _create_synthetic_population(self):
         self.pop_syn_data["person"] = (
             self.pop_syn.loc[:, self.pop_syn_geo_id_columns].join(
                 self.person_stacked_sample))
-        print "\tSize of the housing population table:", (
-            self.pop_syn_data["housing"].shape)
-        print "\tSize of the person population table:", (
-            self.pop_syn_data["person"].shape)
+        print ("\tSize of the housing population table:", (
+            self.pop_syn_data["housing"].shape))
+        print ("\tSize of the person population table:", (
+            self.pop_syn_data["person"].shape))
         # print "Time elapsed for synthetic population 1 is : %.4f" % (
         #    time.time() - t)
 
@@ -214,17 +217,17 @@ def _create_index(self):
         # self.pop_syn_data["housing"].reset_index(inplace=True)
         self.pop_syn_data["housing"].set_index(
             self.pop_syn_housing_matching_id_columns, inplace=True, drop=False)
-        self.pop_syn_data["housing"].sort(inplace=True)
+        self.pop_syn_data["housing"].sort_index(inplace=True)
 
         # self.pop_syn_data["person"].reset_index(inplace=True)
         self.pop_syn_data["person"].set_index(
             self.pop_syn_person_matching_id_columns, inplace=True,
             drop=False)
-        self.pop_syn_data["person"].sort(inplace=True)
+        self.pop_syn_data["person"].sort_index(inplace=True)
         # print "Time elapsed for index is : %.4f" % (time.time() - t)
 
     def export_outputs(self):
-        print "\tGenerating Outputs"
+        print ("\tGenerating Outputs")
         t = time.time()
         self._export_performance_data()
         self._export_multiway_tables()
@@ -232,8 +235,8 @@ def export_outputs(self):
         self._export_weights()
         self._export_synthetic_population()
         self._pretty_print_scenario_configuration_file_to_output()
-        print "\tTime elapsed for generating outputs is : %.4f" % (
-            time.time() - t)
+        print ("\tTime elapsed for generating outputs is : %.4f" % (
+            time.time() - t))
 
     def _pretty_print_scenario_configuration_file_to_output(self):
         filepath = os.path.join(self.outputlocation,
@@ -241,7 +244,11 @@ def _pretty_print_scenario_configuration_file_to_output(self):
         self.scenario_config.write_to_file(filepath)
 
     def _export_performance_data(self):
-        values_to_export = self.scenario_config.outputs.performance
+        try:
+            values_to_export = self.scenario_config.outputs.performance
+        except ConfigError:
+            warn('no performance values to export, set outputs.performance in YAML to generate these')
+            return
         # print "Performance values to export:", values_to_export
         if "ipf" in values_to_export:
             self._export_all_df_in_dict(
@@ -265,7 +272,11 @@ def _export_performance_data(self):
                 self.draw_population_obj.draws_performance, "draws")
 
     def _export_weights(self):
-        export_weights_config = self.scenario_config.outputs.weights
+        try:
+            export_weights_config = self.scenario_config.outputs.weights
+        except ConfigError:
+            warn('no weights to export, set outputs.weights in config.yaml to generate these')
+            return
         if export_weights_config.export:
             df = pd.DataFrame(self.run_ipu_obj.region_sample_weights)
             if export_weights_config.collate_across_geos:
@@ -277,7 +288,7 @@ def _export_df(self, df, filename):
         df.to_csv(filepath)
 
     def _export_all_df_in_dict(self, dict_of_dfs, fileprefix):
-        for key, value in dict_of_dfs.iteritems():
+        for key, value in dict_of_dfs.items():
             filename = "%s%s.csv" % (fileprefix, key)
             filepath = os.path.join(self.outputlocation, filename)
             value.to_csv(filepath)
@@ -285,7 +296,7 @@ def _export_all_df_in_dict(self, dict_of_dfs, fileprefix):
     def _export_multiway_tables(self):
         multiway_tables = self._return_multiway_tables()
 
-        for (filename, filetype), table in multiway_tables.iteritems():
+        for (filename, filetype), table in multiway_tables.items():
             filepath = os.path.join(self.outputlocation, filename)
             table.to_csv(
                 filepath, sep=self.filetype_sep_dict[filetype])
@@ -302,14 +313,18 @@ def _return_multiway_tables(self):
             multiway_table_entity = self._return_aggregate_by_geo(
                 variables, entity_type, entity)
             multiway_tables[(filename, filetype)] = multiway_table_entity
-            print "\t\tTime elapsed for each table is: %.4f" % (
-                time.time() - t)
+            print ("\t\tTime elapsed for each table is: %.4f" % (
+                time.time() - t))
         return multiway_tables
 
     def _export_synthetic_population(self):
         t = time.time()
-        synthetic_population_config = (
-            self.scenario_config.outputs.synthetic_population)
+        try:
+            synthetic_population_config = (
+                self.scenario_config.outputs.synthetic_population)
+        except ConfigError:
+            warn('syn pop files not called, set outputs.synthetic_population in config.yaml to generate these files')
+            return
         sort_columns = self.pop_syn_all_id_columns
         for entity_type in self.entity_types:
             (filename, filetype) = (
@@ -318,14 +333,14 @@ def _export_synthetic_population(self):
             filepath = os.path.join(self.outputlocation, filename)
             # self.pop_syn_data[entity_type].to_csv(
             #    filepath, sep=self.filetype_sep_dict[filetype], index=False)
-            self.pop_syn_data[entity_type].sort(
-                sort_columns, inplace=True)
+            self.pop_syn_data[entity_type].sort_values(
+                by=sort_columns, inplace=True)
             self.pop_syn_data[entity_type].reset_index(drop=True, inplace=True)
             self.pop_syn_data[entity_type].index.name = (
                 "unique_%s_id" % entity_type)
             self.pop_syn_data[entity_type].to_csv(
                 filepath, sep=self.filetype_sep_dict[filetype])
-        print "\tTime to write syn pop files is: %.4f" % (time.time() - t)
+        print ("\tTime to write syn pop files is: %.4f" % (time.time() - t))
 
     def _return_aggregate_by_geo(self, variables, entity_type, entity):
         if isinstance(variables, str):
@@ -335,13 +350,17 @@ def _return_aggregate_by_geo(self, variables, entity_type, entity):
         multiway_table = (self.pop_syn_data[entity_type].groupby(
             groupby_columns).size())
         multiway_table_entity = (
-            multiway_table[entity].unstack(level=range(1, columns_count-1))
+            multiway_table[entity].unstack(level=tuple(range(1, columns_count-1)))
             )
         return multiway_table_entity
 
     def _export_summary(self):
         t = time.time()
-        summary_config = self.scenario_config.outputs.summary
+        try:
+            summary_config = self.scenario_config.outputs.summary
+        except ConfigError:
+            print("\tSummary creation not called in YAML")
+            return
         marginal_geo = self._return_marginal_geo()
         (geo_filename, geo_filetype) = (
             summary_config.geo.filename, summary_config.geo.filetype)
@@ -356,7 +375,7 @@ def _export_summary(self):
         filepath = os.path.join(self.outputlocation, region_filename)
         marginal_region.to_csv(
             filepath, sep=self.filetype_sep_dict[region_filetype])
-        print "\tSummary creation took: %.4f" % (time.time() - t)
+        print ("\tSummary creation took: %.4f" % (time.time() - t))
         # print marginal_region
 
     def _return_marginal_region(self, marginal_geo):
@@ -393,7 +412,7 @@ def _stack_marginal(self, marginal_list):
         stacked_marginal.index.name = "categories"
         stacked_marginal.reset_index(inplace=True)
         stacked_marginal.set_index(["name", "categories"], inplace=True)
-        stacked_marginal.sort(inplace=True)  # Sort by row indices
+        stacked_marginal.sort_index(inplace=True)  # Sort by row indices
         return stacked_marginal.T
 
     def _report_summary(self, geo_id_rows_syn, geo_id_frequencies,
@@ -421,7 +440,7 @@ def _report_summary(self, geo_id_rows_syn, geo_id_frequencies,
         sad_in_constraints = (geo_id_synthetic["abs_diff_constraint"]).sum()
         sd_in_constraints = (geo_id_synthetic["diff_constraint"]).sum()
 
-        print "%.4f, %f, %f, %f, %f, %f" % (stat, p_value, aad_in_frequencies,
+        print ("%.4f, %f, %f, %f, %f, %f" % (stat, p_value, aad_in_frequencies,
                                             aad_in_constraints,
                                             sad_in_constraints,
-                                            sd_in_constraints)
+                                            sd_in_constraints))
diff --git a/popgen/project.py b/popgen/project.py
index c01502c..4f59ed8 100755
--- a/popgen/project.py
+++ b/popgen/project.py
@@ -1,14 +1,15 @@
+from __future__ import print_function
 import logging
 import os
 import time
 import yaml
 
-from config import Config
-from data import DB
-from ipf import Run_IPF
-from reweighting import Run_Reweighting, Reweighting_DS
-from draw import Draw_Population
-from output import Syn_Population
+from .config import Config
+from .data import DB
+from .ipf import Run_IPF
+from .reweighting import Run_Reweighting, Reweighting_DS
+from .draw import Draw_Population
+from .output import Syn_Population
 
 
 class Project(object):
@@ -28,7 +29,8 @@ def _load_config(self):
         # TODO: validating config file for YAML
         # TODO: validating YAML config file for field types
         # TODO: validating YAML for consistency across fields/config elements
-        config_f = file(self.config_loc, "r")
+        print(os.getcwd())
+        config_f = open(self.config_loc, "r")
         config_dict = yaml.load(config_f)
         self._config = Config(config_dict)
         self.column_names_config = self._config.project.inputs.column_names
@@ -47,7 +49,7 @@ def _load_data(self):
     def run_scenarios(self):
         scenarios_config = self._config.project.scenario
         for scenario_config in scenarios_config:
-            print "Running Scenario: %s" % scenario_config.description
+            print ("Running Scenario: %s" % scenario_config.description)
             scenario_obj = Scenario(self.location,
                                     self.entities, self.housing_entities,
                                     self.person_entities,
@@ -84,7 +86,7 @@ def _run_ipf(self):
                                    self.column_names_config,
                                    self.scenario_config, self.db)
         self.run_ipf_obj.run_ipf()
-        print "IPF completed in: %.4f" % (time.time() - self.t)
+        print ("IPF completed in: %.4f" % (time.time() - self.t))
 
     def _run_weighting(self):
         reweighting_config = self.scenario_config.parameters.reweighting
@@ -98,7 +100,7 @@ def _run_weighting(self):
         self.run_reweighting_obj.run_reweighting(
             self.run_ipf_obj.region_constraints,
             self.run_ipf_obj.geo_constraints)
-        print "Reweighting completed in: %.4f" % (time.time() - self.t)
+        print ("Reweighting completed in: %.4f" % (time.time() - self.t))
 
     def _draw_sample(self):
         self.draw_population_obj = Draw_Population(
@@ -109,7 +111,7 @@ def _draw_sample(self):
             self.run_reweighting_obj.geo_stacked,
             self.run_reweighting_obj.region_sample_weights)
         self.draw_population_obj.draw_population()
-        print "Drawing completed in: %.4f" % (time.time() - self.t)
+        print ("Drawing completed in: %.4f" % (time.time() - self.t))
 
     def _report_results(self):
         self.syn_pop_obj = Syn_Population(
@@ -126,11 +128,11 @@ def _report_results(self):
         self.syn_pop_obj.add_records()
         self.syn_pop_obj.prepare_data()
         self.syn_pop_obj.export_outputs()
-        print "Results completed in: %.4f" % (time.time() - self.t)
+        print ("Results completed in: %.4f" % (time.time() - self.t))
 
 
 def popgen_run(project_config):
-    logger = logging.getLgger()
+    logger = logging.getLogger()
     pass
 
 if __name__ == "__main__":
@@ -142,4 +144,4 @@ def popgen_run(project_config):
     p_obj = Project("../tutorials/1_basic_popgen_setup/configuration.yaml")
     p_obj.load_project()
     p_obj.run_scenarios()
-    print "Time it took: %.4f" % (time.time() - t)
+    print ("Time it took: %.4f" % (time.time() - t))
diff --git a/popgen/reweighting.py b/popgen/reweighting.py
old mode 100644
new mode 100755
index e46cd66..3b47ac8
--- a/popgen/reweighting.py
+++ b/popgen/reweighting.py
@@ -1,3 +1,4 @@
+from __future__ import print_function
 import numpy as np
 import pandas as pd
 
@@ -20,7 +21,7 @@ def get_sample_restructure(self, entity, sample, variable_names, hid_name):
         columns_count = len(groupby_columns)
         sample_restruct = (sample.groupby(groupby_columns)
                            .size()
-                           .unstack(level=range(1, columns_count))
+                           .unstack(level=tuple(range(1, columns_count)))
                            .fillna(0)
                            )
         return sample_restruct
@@ -52,11 +53,11 @@ def get_stacked_sample_restruct(self, sample_restruct_list):
             else:
                 stacked_sample = sample_restruct.join(stacked_sample,
                                                       how="outer").fillna(0)
-        stacked_sample.sort(inplace=True)  # Sort by row indices
-        stacked_sample.sort_index(axis=1,
-                                  inplace=True)  # Sort columns alphabetically
+        stacked_sample.sort_index(inplace=True)  # Sort by row indices
         stacked_sample.columns = pd.Index(stacked_sample.columns,
-                                          tuplelize_cols=False)
+                                          tuplelize_cols=True)
+        # stacked_sample.sort_index(axis=1,
+                                  # inplace=True)  # Sort columns alphabetically
         return stacked_sample
 
 
@@ -117,8 +118,8 @@ def _create_reweighting_performance_df(self):
         # TODO: In the future change the frequency at which
         # performance measures are stored as a parameter that is
         # specified by the user
-        self.iters_to_archive = range(0, self.outer_iterations,
-                                      self.archive_performance_frequency)
+        self.iters_to_archive = tuple(range(0, self.outer_iterations,
+                                      self.archive_performance_frequency))
         self.average_diffs = pd.DataFrame(index=self.db.geo_ids,
                                           columns=self.iters_to_archive)
 
@@ -161,7 +162,7 @@ def run_reweighting(self, region_constraints, geo_constraints):
                 #       (time.time() - t))
             self._populate_sample_weights(sample_weights, region_id, geo_ids)
             # print self.average_deviations
-            print "\tsample_weights sum:", sample_weights.sum()
+            print ("\tsample_weights sum:", sample_weights.sum())
 
     def _adjust_sample_weights(self, sample_weights, constraints,
                                iters=1, geo=False):
@@ -236,7 +237,7 @@ def _entropy_adjust_sample_weights(self, sample_weights, constraints,
 
     def _find_equation(self, contrib, weights_mul_contrib):
         root_power_weight = np.bincount(contrib, weights=weights_mul_contrib)
-        root_power = np.array(range(contrib.max() + 1))
+        root_power = np.arange(contrib.max() + 1)
         return root_power[1:], root_power_weight[1:]
 
     def _optimizing_function(self, root, root_power, root_power_weight,
diff --git a/test/test_popgen.py b/test/test_popgen.py
index cafd948..18b0b4b 100755
--- a/test/test_popgen.py
+++ b/test/test_popgen.py
@@ -2,6 +2,12 @@
 Tests for `popgen` module.
 """
 import pytest
+
+import sys
+print('<sys.path>')
+print(sys.path)
+print('</sys.path>')
+
 from popgen import Project