diff --git a/.gitignore b/.gitignore
index eee19081..ccc4a583 100644
--- a/.gitignore
+++ b/.gitignore
@@ -63,7 +63,7 @@ target/
 
 # pyenv python configuration file
 .python-version
-
+.venv
 junit.xml
 
 # crapple
diff --git a/docs/conf.py b/docs/conf.py
index a23474ec..7f8dd108 100644
--- a/docs/conf.py
+++ b/docs/conf.py
@@ -80,7 +80,7 @@
 #
 # This is also used if you do content translation via gettext catalogs.
 # Usually you set "language" from the command line for these cases.
-language = None
+language = 'en'
 
 # List of patterns, relative to source directory, that match files and
 # directories to ignore when looking for source files.
diff --git a/gerrychain/accept.py b/gerrychain/accept.py
index 0884ddc6..c4b02809 100644
--- a/gerrychain/accept.py
+++ b/gerrychain/accept.py
@@ -1,4 +1,13 @@
-from .random import random
+"""
+This module provides the main acceptance function used in ReCom Markov chains.
+
+Dependencies:
+- random: For random number generation for probabilistic acceptance.
+
+Last Updated: 11 Jan 2024
+"""
+
+import random
 from gerrychain.partition import Partition
 
 
@@ -7,12 +16,15 @@ def always_accept(partition: Partition) -> bool:
 
 
 def cut_edge_accept(partition: Partition) -> bool:
-    """Always accepts the flip if the number of cut_edges increases.
+    """
+    Always accepts the flip if the number of cut_edges increases.
     Otherwise, uses the Metropolis criterion to decide.
 
     :param partition: The current partition to accept a flip from.
-    :return: True if accepted, False to remain in place
+    :type partition: Partition
 
+    :return: True if accepted, False to remain in place
+    :rtype: bool
     """
     bound = 1.0
 
diff --git a/gerrychain/chain.py b/gerrychain/chain.py
index 4cf1c34b..8901d282 100644
--- a/gerrychain/chain.py
+++ b/gerrychain/chain.py
@@ -1,3 +1,26 @@
+"""
+This module provides the MarkovChain class, which is designed to facilitate the creation
+and iteration of Markov chains in the context of political redistricting and gerrymandering
+analysis. It allows for the exploration of different districting plans based on specified
+constraints and acceptance criteria.
+
+Key Components:
+- MarkovChain: The main class used for creating and iterating over Markov chain states.
+- Validator: A helper class for validating proposed states in the Markov chain.
+    see :class:`~gerrychain.constraints.Validator` for more details.
+
+Usage:
+The primary use of this module is to create an instance of MarkovChain with appropriate
+parameters like proposal function, constraints, acceptance function, and initial state,
+and then to iterate through the states of the Markov chain, yielding a new proposal
+at each step.
+
+Dependencies:
+- typing: Used for type hints.
+
+Last Updated: 11 Jan 2024
+"""
+
 from .constraints import Validator
 from typing import Union, Iterable, Callable, Optional
 
@@ -7,8 +30,11 @@
 
 class MarkovChain:
     """
-    MarkovChain is an iterator that allows the user to iterate over the states
-    of a Markov chain run.
+    MarkovChain is a class that creates an iterator for iterating over the states
+    of a Markov chain run in a gerrymandering analysis context.
+
+    It allows for the generation of a sequence of partitions (states) of a political
+    districting plans, where each partition represents a possible state in the Markov chain.
 
     Example usage:
 
@@ -30,15 +56,23 @@ def __init__(
     ) -> None:
         """
         :param proposal: Function proposing the next state from the current state.
+        :type proposal: Callable
         :param constraints: A function with signature ``Partition -> bool`` determining whether
             the proposed next state is valid (passes all binary constraints). Usually
             this is a :class:`~gerrychain.constraints.Validator` class instance.
+        :type constraints: Union[Iterable[Callable], Validator, Iterable[Bounds], Callable]
         :param accept: Function accepting or rejecting the proposed state. In the most basic
             use case, this always returns ``True``. But if the user wanted to use a
             Metropolis-Hastings acceptance rule, this is where you would implement it.
+        :type accept: Callable
         :param initial_state: Initial :class:`gerrychain.partition.Partition` class.
+        :type initial_state: Optional[Partition]
         :param total_steps: Number of steps to run.
+        :type total_steps: int
 
+        :return: None
+
+        :raises ValueError: If the initial_state is not valid according to the constraints.
         """
         if callable(constraints):
             is_valid = Validator([constraints])
@@ -65,11 +99,34 @@ def __init__(
         self.state = initial_state
 
     def __iter__(self) -> 'MarkovChain':
+        """
+        Resets the Markov chain iterator.
+
+        This method is called when an iterator is required for a container. It sets the
+        counter to 0 and resets the state to the initial state.
+
+        :return: Returns itself as an iterator object.
+        :rtype: MarkovChain
+        """
         self.counter = 0
         self.state = self.initial_state
         return self
 
     def __next__(self) -> Optional[Partition]:
+        """
+        Advances the Markov chain to the next state.
+
+        This method is called to get the next item in the iteration.
+        It proposes the next state and moves to t it if that state is
+        valid according to the constraints and if accepted by the
+        acceptance function. If the total number of steps has been
+        reached, it raises a StopIteration exception.
+
+        :return: The next state of the Markov chain.
+        :rtype: Optional[Partition]
+
+        :raises StopIteration: If the total number of steps has been reached.
+        """
         if self.counter == 0:
             self.counter += 1
             return self.state
@@ -88,12 +145,26 @@ def __next__(self) -> Optional[Partition]:
         raise StopIteration
 
     def __len__(self) -> int:
+        """
+        Returns the total number of steps in the Markov chain.
+
+        :return: The total number of steps in the Markov chain.
+        :rtype: int
+        """
         return self.total_steps
 
     def __repr__(self) -> str:
         return "<MarkovChain [{} steps]>".format(len(self))
 
     def with_progress_bar(self):
+        """
+        Wraps the Markov chain in a tqdm progress bar.
+
+        Useful for long-running Markov chains where you want to keep track
+        of the progress. Requires the `tqdm` package to be installed.
+
+        :return: A tqdm-wrapped Markov chain.
+        """
         from tqdm.auto import tqdm
 
         return tqdm(self)
diff --git a/gerrychain/constraints/contiguity.py b/gerrychain/constraints/contiguity.py
index 1ec2fa3b..3f58ef16 100644
--- a/gerrychain/constraints/contiguity.py
+++ b/gerrychain/constraints/contiguity.py
@@ -3,7 +3,7 @@
 
 import networkx as nx
 
-from ..random import random
+import random
 from .bounds import SelfConfiguringLowerBound
 
 
diff --git a/gerrychain/grid.py b/gerrychain/grid.py
index 43c34564..b89a3534 100644
--- a/gerrychain/grid.py
+++ b/gerrychain/grid.py
@@ -1,7 +1,18 @@
-import math
+"""
+This module provides a Grid class used for creating and manipulating grid partitions.
+It's part of the GerryChain suite, designed to facilitate experiments with redistricting
+plans without the need for extensive data processing. This module relies on NetworkX for
+graph operations and integrates with GerryChain's Partition class.
+
+Dependencies:
+- math: For math.floor() function.
+- networkx: For graph operations with using the graph structure in
+    :class:`~gerrychain.graph.Graph`.
+- typing: Used for type hints.
+"""
 
+import math
 import networkx
-
 from gerrychain.partition import Partition
 from gerrychain.graph import Graph
 from gerrychain.updaters import (
@@ -14,8 +25,7 @@
     perimeter,
 )
 from gerrychain.metrics import polsby_popper
-
-from typing import Callable, Dict, Optional, Tuple
+from typing import Callable, Dict, Optional, Tuple, Any
 
 
 class Grid(Partition):
@@ -54,15 +64,35 @@ def __init__(
         flips: Optional[Dict[Tuple[int, int], int]] = None,
     ) -> None:
         """
-        :param dimensions: tuple (m,n) of the desired dimensions of the grid.
-        :param with_diagonals: (optional, defaults to False) whether to include diagonals
-            as edges of the graph (i.e., whether to use 'queen' adjacency rather than
-            'rook' adjacency).
-        :param assignment: (optional) dict matching nodes to their districts. If not
-            provided, partitions the grid into 4 quarters of roughly equal size.
-        :param updaters: (optional) dict matching names of attributes of the Partition
-            to functions that compute their values. If not provided, the Grid
-            configures the cut_edges updater for convenience.
+        If the updaters are not specified, the default updaters are used, which are as follows::
+
+            default_updaters = {
+                "cut_edges": cut_edges,
+                "population": Tally("population"),
+                "perimeter": perimeter,
+                "exterior_boundaries": exterior_boundaries,
+                "interior_boundaries": interior_boundaries,
+                "boundary_nodes": boundary_nodes,
+                "area": Tally("area", alias="area"),
+                "polsby_popper": polsby_popper,
+                "cut_edges_by_part": cut_edges_by_part,
+            }
+
+
+        :param dimensions: The grid dimensions (rows, columns), defaults to None.
+        :type dimensions: Tuple[int, int], optional
+        :param with_diagonals: If True, includes diagonal connections, defaults to False.
+        :type with_diagonals: bool, optional
+        :param assignment: Node-to-district assignments, defaults to None.
+        :type assignment: Dict, optional
+        :param updaters: Custom updater functions, defaults to None.
+        :type updaters: Dict[str, Callable], optional
+        :param parent: Parent Grid object for inheritance, defaults to None.
+        :type parent: Grid, optional
+        :param flips: Node flips for partition changes, defaults to None.
+        :type flips: Dict[Tuple[int, int], int], optional
+
+        :raises Exception: If neither dimensions nor parent is provided.
         """
         if dimensions:
             self.dimensions = dimensions
@@ -100,12 +130,32 @@ def as_list_of_lists(self):
         Returns the grid as a list of lists (like a matrix), where the (i,j)th
         entry is the assigned district of the node in position (i,j) on the
         grid.
+
+        :return: List of lists representing the grid.
+        :rtype: List[List[int]]
         """
         m, n = self.dimensions
         return [[self.assignment.mapping[(i, j)] for i in range(m)] for j in range(n)]
 
 
-def create_grid_graph(dimensions: Tuple[int, int], with_diagonals: bool) -> Graph:
+def create_grid_graph(
+    dimensions: Tuple[int, int],
+    with_diagonals: bool
+) -> Graph:
+    """
+    Creates a grid graph with the specified dimensions.
+    Optionally includes diagonal connections between nodes.
+
+    :param dimensions: The grid dimensions (rows, columns).
+    :type dimensions: Tuple[int, int]
+    :param with_diagonals: If True, includes diagonal connections.
+    :type with_diagonals: bool
+
+    :return: A grid graph.
+    :rtype: Graph
+
+    :raises ValueError: If the dimensions are not a tuple of length 2.
+    """
     if len(dimensions) != 2:
         raise ValueError("Expected two dimensions.")
     m, n = dimensions
@@ -133,12 +183,43 @@ def create_grid_graph(dimensions: Tuple[int, int], with_diagonals: bool) -> Grap
     return graph
 
 
-def give_constant_attribute(graph, attribute, value):
+def give_constant_attribute(
+    graph: Graph,
+    attribute: Any,
+    value: Any
+) -> None:
+    """
+    Sets the specified attribute to the specified value for all nodes in the graph.
+
+    :param graph: The graph to modify.
+    :type graph: Graph
+    :param attribute: The attribute to set.
+    :type attribute: Any
+    :param value: The value to set the attribute to.
+    :type value: Any
+
+    :return: None
+    """
     for node in graph.nodes:
         graph.nodes[node][attribute] = value
 
 
-def tag_boundary_nodes(graph: Graph, dimensions: Tuple[int, int]) -> None:
+def tag_boundary_nodes(
+    graph: Graph,
+    dimensions: Tuple[int, int]
+) -> None:
+    """
+    Adds the boolean attribute ``boundary_node`` to each node in the graph.
+    If the node is on the boundary of the grid, that node also gets the attribute
+    ``boundary_perim`` which is determined by the function :func:`get_boundary_perim`.
+
+    :param graph: The graph to modify.
+    :type graph: Graph
+    :param dimensions: The dimensions of the grid.
+    :type dimensions: Tuple[int, int]
+
+    :return: None
+    """
     m, n = dimensions
     for node in graph.nodes:
         if node[0] in [0, m - 1] or node[1] in [0, n - 1]:
@@ -148,7 +229,23 @@ def tag_boundary_nodes(graph: Graph, dimensions: Tuple[int, int]) -> None:
             graph.nodes[node]["boundary_node"] = False
 
 
-def get_boundary_perim(node: Tuple[int, int], dimensions: Tuple[int, int]) -> int:
+def get_boundary_perim(
+    node: Tuple[int, int],
+    dimensions: Tuple[int, int]
+) -> int:
+    """
+    Determines the boundary perimeter of a node on the grid.
+    The boundary perimeter is the number of sides of the node that
+    are on the boundary of the grid.
+
+    :param node: The node to check.
+    :type node: Tuple[int, int]
+    :param dimensions: The dimensions of the grid.
+    :type dimensions: Tuple[int, int]
+
+    :return: The boundary perimeter of the node.
+    :rtype: int
+    """
     m, n = dimensions
     if node in [(0, 0), (m - 1, 0), (0, n - 1), (m - 1, n - 1)]:
         return 2
@@ -158,7 +255,7 @@ def get_boundary_perim(node: Tuple[int, int], dimensions: Tuple[int, int]) -> in
         return 0
 
 
-def color_half(node, threshold=5):
+def color_half(node: Tuple[int, int], threshold: int = 5) -> int:
     x = node[0]
     return 0 if x <= threshold else 1
 
@@ -168,19 +265,3 @@ def color_quadrants(node: Tuple[int, int], thresholds: Tuple[int, int]) -> int:
     x_color = 0 if x < thresholds[0] else 1
     y_color = 0 if y < thresholds[1] else 2
     return x_color + y_color
-
-
-def grid_size(parition):
-    """ This is a hardcoded population function
-    for the grid class"""
-
-    L = parition.as_list_of_lists()
-    permit = [3, 4, 5]
-
-    sizes = [0, 0, 0, 0]
-
-    for i in range(len(L)):
-        for j in range(len(L[0])):
-            sizes[L[i][j]] += 1
-
-    return all(x in permit for x in sizes)
diff --git a/gerrychain/metagraph.py b/gerrychain/metagraph.py
index ffbc28bc..6a02fbb3 100644
--- a/gerrychain/metagraph.py
+++ b/gerrychain/metagraph.py
@@ -1,12 +1,31 @@
-from itertools import product
+"""
+This module provides the main tools for interacting with the metagraph of partitions.
+The metagraph of partitions is the set of partitions that are reachable from the
+current partition by a single flip.
+
+Dependencies:
+- itertools: Used for product() function.
+- typing: Used for type hints.
 
+Last Updated: 11 Jan 2024
+"""
+
+from itertools import product
 from .constraints import Validator
 from typing import Callable, Dict, Iterator, Iterable, Union
-
 from gerrychain.partition import Partition
 
 
 def all_cut_edge_flips(partition: Partition) -> Iterator[Dict]:
+    """
+    Generate all possible flips of cut edges in a partition
+    without any contraints.
+
+    :param partition: The partition object.
+    :type partition: Partition
+    :return: An iterator that yields dictionaries representing the flipped edges.
+    :rtype: Iterator[Dict]
+    """
     for edge, index in product(partition.cut_edges, (0, 1)):
         yield {edge[index]: partition.assignment.mapping[edge[1 - index]]}
 
@@ -15,9 +34,21 @@ def all_valid_states_one_flip_away(
     partition: Partition,
     constraints: Union[Iterable[Callable], Callable]
 ) -> Iterator[Partition]:
-    """Generates all valid Partitions that differ from the given partition
+    """
+    Generates all valid Partitions that differ from the given partition
     by one flip. These are the given partition's neighbors in the metagraph
-    of partitions.
+    of partitions. (The metagraph of partitions is the set of partitions
+    that is reachable from the given partition by a single flip under the
+    prescribed constraints.)
+
+    :param partition: The initial partition.
+    :type partition: Partition
+    :param constraints: Constraints to determine the validity of a partition.
+                        It can be a single callable or an iterable of callables.
+    :type constraints: Union[Iterable[Callable], Callable]
+    :return: An iterator that yields all valid partitions that differ from the
+             given partition by one flip.
+    :rtype: Iterator[Partition]
     """
     if callable(constraints):
         is_valid = constraints
@@ -34,6 +65,17 @@ def all_valid_flips(
     partition: Partition,
     constraints: Union[Iterable[Callable], Callable]
 ) -> Iterator[Dict]:
+    """
+    Generate all valid flips for a given partition subject
+    to the prescribed constraints.
+
+    :param partition: The initial partition.
+    :type partition: Partition
+    :param constraints: The constraints to be satisfied.
+    :type constraints: Union[Iterable[Callable], Callable]
+    :return: An iterator that yields dictionaries representing valid flips.
+    :rtype: Iterator[Dict]
+    """
     for state in all_valid_states_one_flip_away(partition, constraints):
         yield state.flips
 
@@ -42,4 +84,18 @@ def metagraph_degree(
     partition: Partition,
     constraints: Union[Iterable[Callable], Callable]
 ) -> int:
+    """
+    Calculate the degree of the metagraph for a given partition.
+    That is to say, compute how many possible valid states are reachable from
+    the state given by partition in a single flip subject to the prescribed
+    constraints.
+
+    :param partition: The partition object representing the current state.
+    :type partition: Partition
+    :param constraints: The constraints to be applied to the partition.
+                        It can be a single constraint or an iterable of constraints.
+    :type constraints: Union[Iterable[Callable], Callable]
+    :return: The degree of the metagraph.
+    :rtype: int
+    """
     return len(list(all_valid_states_one_flip_away(partition, constraints)))
diff --git a/gerrychain/proposals/proposals.py b/gerrychain/proposals/proposals.py
index b120d919..975c95af 100644
--- a/gerrychain/proposals/proposals.py
+++ b/gerrychain/proposals/proposals.py
@@ -1,4 +1,4 @@
-from ..random import random
+import random
 
 
 def propose_any_node_flip(partition):
diff --git a/gerrychain/proposals/spectral_proposals.py b/gerrychain/proposals/spectral_proposals.py
index 09048955..77664dcf 100644
--- a/gerrychain/proposals/spectral_proposals.py
+++ b/gerrychain/proposals/spectral_proposals.py
@@ -1,6 +1,6 @@
 import networkx as nx
 from numpy import linalg as LA
-from ..random import random
+import random
 
 
 def spectral_cut(graph, part_labels, weight_type, lap_type):
diff --git a/gerrychain/proposals/tree_proposals.py b/gerrychain/proposals/tree_proposals.py
index 2d2cb684..4bb781f8 100644
--- a/gerrychain/proposals/tree_proposals.py
+++ b/gerrychain/proposals/tree_proposals.py
@@ -1,29 +1,29 @@
 from functools import partial
 from inspect import signature
-from ..random import random
+import random
 
+from gerrychain.partition import Partition
 from ..tree import (
     recursive_tree_part, bipartition_tree, bipartition_tree_random,
     _bipartition_tree_random_all, uniform_spanning_tree,
     find_balanced_edge_cuts_memoization,
 )
+from typing import Callable, Optional, Dict
 
 
 def recom(
-    partition, pop_col, pop_target, epsilon, node_repeats=1, 
-    weight_dict = None,
-    method=bipartition_tree
-):
-    """ReCom proposal.
-
-    Description from MGGG's 2018 Virginia House of Delegates report:
-    At each step, we (uniformly) randomly select a pair of adjacent districts and
-    merge all of their blocks in to a single unit. Then, we generate a spanning tree
-    for the blocks of the merged unit with the Kruskal/Karger algorithm. Finally,
-    we cut an edge of the tree at random, checking that this separates the region
-    into two new districts that are population balanced.
+    partition: Partition,
+    pop_col: str,
+    pop_target: float,
+    epsilon: float,
+    node_repeats: int = 1,
+    weight_dict: Optional[Dict] = None,
+    method: Callable = bipartition_tree
+) -> Partition:
+    """
+    Example usage:
 
-    Example usage::
+    .. code-block:: python
 
         from functools import partial
         from gerrychain import MarkovChain
@@ -41,6 +41,7 @@ def recom(
         chain = MarkovChain(proposal, constraints, accept, partition, total_steps)
 
     """
+
     edge = random.choice(tuple(partition["cut_edges"]))
     parts_to_merge = (partition.assignment.mapping[edge[0]], partition.assignment.mapping[edge[1]])
 
diff --git a/gerrychain/random.py b/gerrychain/random.py
deleted file mode 100644
index 2a1940e1..00000000
--- a/gerrychain/random.py
+++ /dev/null
@@ -1,5 +0,0 @@
-import os
-import random
-
-seed = os.environ.get("GERRYCHAIN_RANDOM_SEED", 2018)
-random.seed(seed)
diff --git a/gerrychain/tree.py b/gerrychain/tree.py
index 1ccc2d8a..4924a7ff 100644
--- a/gerrychain/tree.py
+++ b/gerrychain/tree.py
@@ -1,9 +1,38 @@
+"""
+This module provides tools and algorithms for manipulating and analyzing graphs,
+particularly focused on partitioning graphs based on population data. It leverages the
+NetworkX library to handle graph structures and implements various algorithms for graph
+partitioning and tree traversal.
+
+Key functionalities include:
+
+- Predecessor and successor functions for graph traversal using breadth-first search.
+- Implementation of random and uniform spanning trees for graph partitioning.
+- The `PopulatedGraph` class, which represents a graph with additional population data,
+  and methods for assessing and modifying this data.
+- Functions for finding balanced edge cuts in a populated graph, either through
+  contraction or memoization techniques.
+- A suite of functions (`bipartition_tree`, `recursive_tree_part`, `get_seed_chunks`, etc.)
+  for partitioning graphs into balanced subsets based on population targets and tolerances.
+- Utility functions like `get_max_prime_factor_less_than` and `recursive_seed_part_inner`
+  to assist in complex partitioning tasks.
+
+Dependencies:
+
+- networkx: Used for graph data structure and algorithms.
+- random: Provides random number generation for probabilistic approaches.
+- typing: Used for type hints.
+
+Last Updated: 11 Jan 2024
+"""
+
+
 import networkx as nx
 from networkx.algorithms import tree
 
 from functools import partial
 from inspect import signature
-from .random import random
+import random
 from collections import deque, namedtuple
 from typing import Any, Callable, Dict, List, Optional, Set, Union, Hashable, Sequence, Tuple
 
@@ -16,27 +45,28 @@ def successors(h: nx.Graph, root: Any) -> Dict:
     return {a: b for a, b in nx.bfs_successors(h, root)}
 
 
-def random_spanning_tree(graph: nx.Graph, weight_dict: Dict) -> nx.Graph:
-    """ 
+def random_spanning_tree(graph: nx.Graph, weight_dict: Optional[Dict] = None) -> nx.Graph:
+    """
     Builds a spanning tree chosen by Kruskal's method using random weights.
-    
+
     :param graph: The input graph to build the spanning tree from. Should be a Networkx Graph.
     :type graph: nx.Graph
-    :param weight_dict: Dictionary of weights to add to the random weights used in region-aware variants.
-    :type weight_dict: Dict
+    :param weight_dict: Dictionary of weights to add to the random weights used in region-aware
+        variants.
+    :type weight_dict: Optional[Dict], optional
     :return: The maximal spanning tree represented as a Networkx Graph.
     :rtype: nx.Graph
     """
     if weight_dict is None:
         weight_dict = dict()
-        
+
     for edge in graph.edges():
         weight = random.random()
         for key, value in weight_dict.items():
             if graph.nodes[edge[0]][key] == graph.nodes[edge[1]][key] and \
                graph.nodes[edge[0]][key] is not None:
                 weight += value
-                
+
         graph.edges[edge]["random_weight"] = weight
 
     spanning_tree = tree.maximum_spanning_tree(
@@ -45,13 +75,14 @@ def random_spanning_tree(graph: nx.Graph, weight_dict: Dict) -> nx.Graph:
     return spanning_tree
 
 
-def uniform_spanning_tree( 
-    graph: nx.Graph, 
+def uniform_spanning_tree(
+    graph: nx.Graph,
     choice: Callable = random.choice
 ) -> nx.Graph:
-    """ 
+    """
     Builds a spanning tree chosen uniformly from the space of all
     spanning trees of the graph. Uses Wilson's algorithm.
+
     :param graph: Networkx Graph
     :type graph: nx.Graph
     :param choice: :func:`random.choice`. Defaults to :func:`random.choice`.
@@ -90,10 +121,11 @@ class PopulatedGraph:
     :type populations: Dict
     :param ideal_pop: The ideal population for each district.
     :type ideal_pop: float
-    :param epsilon: The tolerance for population deviation from the ideal population within each 
-    district.
+    :param epsilon: The tolerance for population deviation from the ideal population within each
+        district.
     :type epsilon: float
     """
+
     def __init__(
         self,
         graph: nx.Graph,
@@ -125,15 +157,26 @@ def has_ideal_population(self, node) -> bool:
             abs(self.population[node] - self.ideal_pop) < self.epsilon * self.ideal_pop
         )
 
+    def __repr__(self) -> str:
+        graph_info = f"Graph(nodes={len(self.graph.nodes)}, edges={len(self.graph.edges)})"
+        return (
+            f"{self.__class__.__name__}("
+            f"graph={graph_info}, "
+            f"total_population={self.tot_pop}, "
+            f"ideal_population={self.ideal_pop}, "
+            f"epsilon={self.epsilon})"
+        )
 
 
 # Tuple that is used in the find_balanced_edge_cuts function
-# Comment added to make this easier to find
 Cut = namedtuple("Cut", "edge subset")
+Cut.__doc__ = "Represents a cut in a graph."
+Cut.edge.__doc__ = "The edge where the cut is made."
+Cut.subset.__doc__ = "The subset of nodes on one side of the cut."
 
 
 def find_balanced_edge_cuts_contraction(
-    h: PopulatedGraph, 
+    h: PopulatedGraph,
     choice: Callable = random.choice
 ) -> List[Cut]:
     """
@@ -172,9 +215,10 @@ def find_balanced_edge_cuts_memoization(
     """
     Find balanced edge cuts using memoization.
 
-    This function takes a PopulatedGraph object and a choice function as input and returns a list of balanced edge cuts.
-    A balanced edge cut is defined as a cut that divides the graph into two subsets, such that the population of each subset
-    is close to the ideal population defined by the PopulatedGraph object.
+    This function takes a PopulatedGraph object and a choice function as input and returns a list
+    of balanced edge cuts. A balanced edge cut is defined as a cut that divides the graph into
+    two subsets, such that the population of each subset is close to the ideal population
+    defined by the PopulatedGraph object.
 
     :param h: The PopulatedGraph object representing the graph.
     :type h: PopulatedGraph
@@ -183,7 +227,7 @@ def find_balanced_edge_cuts_memoization(
     :return: A list of balanced edge cuts.
     :rtype: List[Any]
     """
-    
+
     root = choice([x for x in h if h.degree(x) > 1])
     pred = predecessors(h.graph, root)
     succ = successors(h.graph, root)
@@ -238,7 +282,7 @@ def bipartition_tree(
     node_repeats: int = 1,
     spanning_tree: Optional[nx.Graph] = None,
     spanning_tree_fn: Callable = random_spanning_tree,
-    weight_dict: Dict = None,
+    weight_dict: Optional[Dict] = None,
     balance_edge_fn: Callable = find_balanced_edge_cuts_memoization,
     choice: Callable = random.choice,
     max_attempts: Optional[int] = 10000
@@ -252,38 +296,37 @@ def bipartition_tree(
     Builds up a connected subgraph with a connected complement whose population
     is ``epsilon * pop_target`` away from ``pop_target``.
 
-
     :param graph: The graph to partition.
     :type graph: nx.Graph
     :param pop_col: The node attribute holding the population of each node.
     :type pop_col: str
     :param pop_target: The target population for the returned subset of nodes.
     :type pop_target: Union[int, float]
-    :param epsilon: The allowable deviation from ``pop_target`` (as a percentage of 
+    :param epsilon: The allowable deviation from ``pop_target`` (as a percentage of
         ``pop_target``) for the subgraph's population.
     :type epsilon: float
-    :param node_repeats: A parameter for the algorithm: how many different choices 
+    :param node_repeats: A parameter for the algorithm: how many different choices
         of root to use before drawing a new spanning tree. Defaults to 1.
     :type node_repeats: int
-    :param spanning_tree: The spanning tree for the algorithm to use (used when the 
+    :param spanning_tree: The spanning tree for the algorithm to use (used when the
         algorithm chooses a new root and for testing).
     :type spanning_tree: Optional[nx.Graph]
-    :param spanning_tree_fn: The random spanning tree algorithm to use if a spanning 
+    :param spanning_tree_fn: The random spanning tree algorithm to use if a spanning
         tree is not provided. Defaults to :func:`random_spanning_tree`.
     :type spanning_tree_fn: Callable
-    :param weight_dict: A dictionary of weights for the spanning tree algorithm. 
+    :param weight_dict: A dictionary of weights for the spanning tree algorithm.
         Defaults to None.
-    :type weight_dict: Dict, optional
-    :param balance_edge_fn: The function to find balanced edge cuts. Defaults to 
+    :type weight_dict: Optional[Dict], optional
+    :param balance_edge_fn: The function to find balanced edge cuts. Defaults to
         :func:`find_balanced_edge_cuts_memoization`.
     :type balance_edge_fn: Callable, optional
     :param choice: The function to make a random choice. Can be substituted for testing.
         Defaults to :func:`random.choice`.
     :type choice: Callable
-    :param max_attempts: The maximum number of attempts that should be made to bipartition. 
+    :param max_attempts: The maximum number of attempts that should be made to bipartition.
         Defaults to 1000.
     :type max_attempts: Optional[int]
-    :return: A subset of nodes of ``graph`` (whose induced subgraph is connected). The other 
+    :return: A subset of nodes of ``graph`` (whose induced subgraph is connected). The other
         part of the partition is the complement of this subset.
     :rtype: Set
     :raises RuntimeError: If a possible cut cannot be found after the maximum number of attempts.
@@ -291,7 +334,7 @@ def bipartition_tree(
     # Try to add the region-aware in if the spanning_tree_fn accepts a weight dictionary
     if 'weight_dict' in signature(spanning_tree_fn).parameters:
         spanning_tree_fn = partial(spanning_tree_fn, weight_dict=weight_dict)
-    
+
     populations = {node: graph.nodes[node][pop_col] for node in graph.node_indices}
 
     possible_cuts = []
@@ -342,26 +385,31 @@ def _bipartition_tree_random_all(
     :type epsilon: float
     :param node_repeats: The number of times to repeat the bipartitioning process. Defaults to 1.
     :type node_repeats: int, optional
-    :param repeat_until_valid: Whether to repeat the bipartitioning process until a valid bipartition is found. Defaults to True.
+    :param repeat_until_valid: Whether to repeat the bipartitioning process until a valid
+        bipartition is found. Defaults to True.
     :type repeat_until_valid: bool, optional
-    :param spanning_tree: The spanning tree to use for bipartitioning. If None, a random spanning tree will be generated. Defaults to None.
+    :param spanning_tree: The spanning tree to use for bipartitioning. If None, a random spanning
+        tree will be generated. Defaults to None.
     :type spanning_tree: Optional[nx.Graph], optional
-    :param spanning_tree_fn: The function to generate a spanning tree. Defaults to random_spanning_tree.
+    :param spanning_tree_fn: The function to generate a spanning tree. Defaults to
+        random_spanning_tree.
     :type spanning_tree_fn: Callable, optional
-    :param balance_edge_fn: The function to find balanced edge cuts. Defaults to find_balanced_edge_cuts_memoization.
+    :param balance_edge_fn: The function to find balanced edge cuts. Defaults to
+        find_balanced_edge_cuts_memoization.
     :type balance_edge_fn: Callable, optional
     :param choice: The function to choose a random element from a list. Defaults to random.choice.
     :type choice: Callable, optional
-    :param max_attempts: The maximum number of attempts to find a valid bipartition. If None, there is no limit. Defaults to None.
+    :param max_attempts: The maximum number of attempts to find a valid bipartition. If None,
+        there is no limit. Defaults to None.
     :type max_attempts: Optional[int], optional
 
     :returns: A list of possible cuts that bipartition the tree into two subgraphs.
     :rtype: List[Tuple[Hashable, Hashable]]
 
-    :raises RuntimeError: If a valid bipartition cannot be found after the specified number of attempts.
+    :raises RuntimeError: If a valid bipartition cannot be found after the specified number of
+        attempts.
     """
 
-    
     populations = {node: graph.nodes[node][pop_col] for node in graph.node_indices}
 
     possible_cuts = []
@@ -414,37 +462,38 @@ def bipartition_tree_random(
     Builds up a connected subgraph with a connected complement whose population
     is ``epsilon * pop_target`` away from ``pop_target``.
 
-    :param graph: The graph to partition (must be an instance of nx.Graph)
+    :param graph: The graph to partition
     :type graph: nx.Graph
-    :param pop_col: The node attribute holding the population of each node (must be a string)
+    :param pop_col: The node attribute holding the population of each node
     :type pop_col: str
-    :param pop_target: The target population for the returned subset of nodes (must be an int or float)
+    :param pop_target: The target population for the returned subset of nodes
     :type pop_target: Union[int, float]
     :param epsilon: The allowable deviation from  ``pop_target`` (as a percentage of
-        ``pop_target``) for the subgraph's population (must be a float)
+        ``pop_target``) for the subgraph's population
     :type epsilon: float
     :param node_repeats: A parameter for the algorithm: how many different choices
-        of root to use before drawing a new spanning tree (default is 1, must be an int)
+        of root to use before drawing a new spanning tree. Defaults to 1.
     :type node_repeats: int
     :param repeat_until_valid: Determines whether to keep drawing spanning trees
         until a tree with a balanced cut is found. If `True`, a set of nodes will
         always be returned; if `False`, `None` will be returned if a valid spanning
-        tree is not found on the first try (default is True, must be a bool)
+        tree is not found on the first try. Defaults to True.
     :type repeat_until_valid: bool
     :param spanning_tree: The spanning tree for the algorithm to use (used when the
-        algorithm chooses a new root and for testing) (must be an instance of nx.Graph or None)
+        algorithm chooses a new root and for testing)
     :type spanning_tree: Optional[nx.Graph]
     :param spanning_tree_fn: The random spanning tree algorithm to use if a spanning
-        tree is not provided (must be a callable)
+        tree is not provided
     :type spanning_tree_fn: Callable
-    :param balance_edge_fn: The algorithm used to find balanced cut edges (must be a callable)
+    :param balance_edge_fn: The algorithm used to find balanced cut edges
     :type balance_edge_fn: Callable
-    :param choice: :func:`random.choice`. Can be substituted for testing. (must be a callable)
+    :param choice: :func:`random.choice`. Can be substituted for testing.
     :type choice: Callable
-    :param max_attempts: The max number of attempts that should be made to bipartition. (must be an int or None)
+    :param max_attempts: The max number of attempts that should be made to bipartition.
     :type max_attempts: Optional[int]
 
-    :return: A subset of nodes of ``graph`` (whose induced subgraph is connected) or None if a valid spanning tree is not found.
+    :return: A subset of nodes of ``graph`` (whose induced subgraph is connected) or None if a
+        valid spanning tree is not found.
     :rtype: Union[Set[Any], None]
     """
     possible_cuts = _bipartition_tree_random_all(
@@ -492,9 +541,10 @@ def recursive_tree_part(
     :param node_repeats: Parameter for :func:`~gerrychain.tree_methods.bipartition_tree` to use.
         Defaluts to 1.
     :type node_repeats: int, optional
-    :param method: The partition method to use. Defaults to 
+    :param method: The partition method to use. Defaults to
         `partial(bipartition_tree, max_attempts=10000)`.
     :type method: Callable, optional
+
     :return: New assignments for the nodes of ``graph``.
     :rtype: dict
     """
@@ -566,12 +616,13 @@ def get_seed_chunks(
     :param epsilon: How far (as a percentage of ``pop_target``) from ``pop_target`` the parts
         of the partition can be
     :type epsilon: float
-    :param node_repeats: Parameter for :func:`~gerrychain.tree_methods.bipartition_tree_random` 
+    :param node_repeats: Parameter for :func:`~gerrychain.tree_methods.bipartition_tree_random`
         to use.
     :type node_repeats: int, optional
     :param method: The method to use for bipartitioning the graph.
         Defaults to :func:`~gerrychain.tree_methods.bipartition_tree_random`
     :type method: Callable, optional
+
     :return: New assignments for the nodes of ``graph``.
     :rtype: dict
     """
@@ -651,14 +702,16 @@ def get_max_prime_factor_less_than(
     n: int, ceil: int
 ) -> Optional[int]:
     """
-    Helper function for recursive_seed_part_inner. Returns the largest prime factor of ``n`` less than
-    ``ceil``, or None if all are greater than ceil.
+    Helper function for recursive_seed_part_inner. Returns the largest prime factor of ``n``
+        less than ``ceil``, or None if all are greater than ceil.
 
     :param n: The number to find the largest prime factor for.
     :type n: int
     :param ceil: The upper limit for the largest prime factor.
     :type ceil: int
-    :return: The largest prime factor of ``n`` less than ``ceil``, or None if all are greater than ceil.
+
+    :return: The largest prime factor of ``n`` less than ``ceil``, or None if all are greater
+        than ceil.
     :rtype: int or None
     """
     if n <= 1 or ceil <= 1:
@@ -668,7 +721,7 @@ def get_max_prime_factor_less_than(
     while n % 2 == 0:
         largest_factor = 2
         n //= 2
-    
+
     i = 3
     while i * i <= n:
         while n % i == 0:
@@ -681,7 +734,6 @@ def get_max_prime_factor_less_than(
         largest_factor = n
 
     return largest_factor
- 
 
 
 def recursive_seed_part_inner(
@@ -702,14 +754,16 @@ def recursive_seed_part_inner(
     Splits graph into num_chunks chunks, and then recursively splits each chunk into
     ``num_dists``/num_chunks chunks.
     The number num_chunks of chunks is chosen based on ``n`` and ``ceil`` as follows:
-        If ``n`` is None, and ``ceil`` is None, num_chunks is the largest prime factor
-        of ``num_dists``.
-        If ``n`` is None and ``ceil`` is an integer at least 2, then num_chunks is the
-        largest prime factor of ``num_dists`` that is less than ``ceil``
-        If ``n`` is a positive integer, num_chunks equals n.
+
+    - If ``n`` is None, and ``ceil`` is None, num_chunks is the largest prime factor
+      of ``num_dists``.
+    - If ``n`` is None and ``ceil`` is an integer at least 2, then num_chunks is the
+      largest prime factor of ``num_dists`` that is less than ``ceil``
+    - If ``n`` is a positive integer, num_chunks equals n.
+
     Finally, if the number of chunks as chosen above does not divide ``num_dists``, then
     this function bites off a single district from the graph and recursively partitions
-    the remaining graph into ``num_dists``-1 districts.
+    the remaining graph into ``num_dists - 1`` districts.
 
     :param graph: The graph
     :param num_dists: number of districts to partition the graph into
@@ -728,9 +782,11 @@ def recursive_seed_part_inner(
         If ``ceil`` is a positive integer then finds the largest factor of ``num_dists`` less
         than or equal to ``ceil``, and recursively splits graph into that number of chunks, or
         bites off a district if that number is 1.
+
     :return: New assignments for the nodes of ``graph``.
     :rtype: List of lists, each list is a district
     """
+
     # Chooses num_chunks
     if n is None:
         if ceil is None:
@@ -760,13 +816,13 @@ def recursive_seed_part_inner(
         )
         remaining_nodes -= nodes
         assignment = [nodes] + recursive_seed_part_inner(graph.subgraph(remaining_nodes),
-            num_dists - 1,
-            pop_target,
-            pop_col,
-            epsilon,
-            method,
-            n=n,
-            ceil=ceil)
+                                                         num_dists - 1,
+                                                         pop_target,
+                                                         pop_col,
+                                                         epsilon,
+                                                         method,
+                                                         n=n,
+                                                         ceil=ceil)
 
     # split graph into num_chunks chunks, and recurse into each chunk
     elif num_dists % num_chunks == 0:
@@ -840,6 +896,7 @@ def recursive_seed_part(
         equal to ``ceil``, and recursively splits graph into that number of chunks, or bites off a
         district if that number is 1. Defaults to None.
     :type ceil: Optional[int]
+
     :return: New assignments for the nodes of ``graph``.
     :rtype: dict
     """
diff --git a/readthedocs.yml b/readthedocs.yml
index 599083ee..2b39a761 100644
--- a/readthedocs.yml
+++ b/readthedocs.yml
@@ -1,6 +1,23 @@
+# .readthedocs.yaml
+# Read the Docs configuration file
+# See https://docs.readthedocs.io/en/stable/config-file/v2.html for details
+
+# Required
+version: 2
+
+# OS version to use on the build environment
 build:
-  image: latest
+  os: ubuntu-22.04
+  tools:
+    python: "3.9"
+
+
+# Build documentation in the docs/ directory with Sphinx
+sphinx:
+   builder: dirhtml
+   configuration: docs/conf.py
 
+# Optionally set the version of Python and requirements required to build your docs
 python:
-  version: 3.9
-  setup_py_install: true
+   install:
+   - requirements: docs/requirements.txt
\ No newline at end of file
diff --git a/tests/conftest.py b/tests/conftest.py
index 96d72dc0..da30c585 100644
--- a/tests/conftest.py
+++ b/tests/conftest.py
@@ -1,9 +1,10 @@
 import pytest
 
 from gerrychain import Graph, Partition
-from gerrychain.random import random
+import random
 from gerrychain.updaters import cut_edges
 import networkx
+random.seed(2018)
 
 
 @pytest.fixture
diff --git a/tests/test_tally.py b/tests/test_tally.py
index ecac2dcf..220bc149 100644
--- a/tests/test_tally.py
+++ b/tests/test_tally.py
@@ -5,9 +5,9 @@
 from gerrychain.constraints import no_vanishing_districts, single_flip_contiguous
 from gerrychain.grid import Grid
 from gerrychain.proposals import propose_random_flip
-from gerrychain.random import random
+import random
 from gerrychain.updaters.tally import DataTally, Tally
-
+random.seed(2018)
 
 def random_assignment(graph, num_districts):
     return {node: random.choice(range(num_districts)) for node in graph.nodes}
diff --git a/tests/updaters/test_updaters.py b/tests/updaters/test_updaters.py
index 20e2136a..37a4b97e 100644
--- a/tests/updaters/test_updaters.py
+++ b/tests/updaters/test_updaters.py
@@ -8,13 +8,13 @@
 from gerrychain.graph import Graph
 from gerrychain.partition import Partition
 from gerrychain.proposals import propose_random_flip
-from gerrychain.random import random
+import random
 from gerrychain.updaters import (Election, Tally, boundary_nodes, cut_edges,
                                  cut_edges_by_part, exterior_boundaries,
                                  exterior_boundaries_as_a_set,
                                  interior_boundaries, perimeter)
 from gerrychain.updaters.election import ElectionResults
-
+random.seed(2018)
 
 @pytest.fixture
 def graph_with_d_and_r_cols(graph_with_random_data_factory):