How to Fix the House of Representatives

The representation in Congress is poorly distributed across states. How can we improve this?

Loading the Data

While Mathematica has population conveniently built in to state entities through the Wolfram Data Repository, let’s be extra careful and load the same decennial Census figures the government has used for apportionment, as well as the official number of representatives per state allotted every tens years. The historical data is drawn from Census.gov in the Notebook linked at the top of this Notebook and published to the Cloud:
(HyperlinkURL["Download data on Congressional apportionment and the decennial Census", https://www.wolframcloud.com/objects/christopher.e.wilson/Published/getPopulationApportionment.nb])

Clear["Global`*"];
censusData = CloudGet["https://www.wolframcloud.com/objects/1a4bfb31-d5d4-49e3-9e04-4cf80e1fae1b"];
stateAbbrs = Keys@censusData;
(* We need a list of state abbreviations without DC since it is not factored into allocation -- "Taxation Without Representation!" *)
stateAbbrs50 = Select[stateAbbrs, # != "DC"&];

How apportionment works: The Huntington-Hill method

The Constitution stipulates that the number of representatives per state should be reapportioned after each decennial census based on state populations. It does not specify the total number of representatives or the exact method of allotment, other than specifying that each state should have at least one representative.

Since the 1940 Census, the United States has used the Huntington-Hill algorithm, known as the "method of equal proportions," to divvy up the 435 House seats. This begins by allotting one seat to each state, leaving 385 remaining. The 50 states are then allocated these remaining seats one-by-one according to which one has the highest priority according to a simple equation: its population divided by the square root of its current seats times one extra seat (the Geometric Mean):

(P/(Sqrt[n (1+n)]))

Lets recreate this method

HuntingtonHill[stateAbbr_, stateReps_] := 1.0 * censusData[stateAbbr]["population_decennial"]["2010"] / Sqrt[stateReps * (stateReps + 1)]

The allocation function will assign electoral votes according to the above formula one at a time. To save computational time, we’ll store the priorities until we update the number of seats in a state rather than run 50 square roots 385 times.

ApportionOne[reps_, priorities_] := (
	topPriority = First[Keys[priorities]];
	updatedPriorities = priorities;
	updatedReps = reps;
	updatedReps[topPriority] += 1;
	updatedPriorities[topPriority] = HuntingtonHill[topPriority, updatedReps[topPriority]];
	{ updatedReps, Reverse@Sort@updatedPriorities }
);

Since we might want to tinker with the size of the House, let’s make that a parameter, total_, in the function that runs the entire apportionment. It will return an association with the final apportionment as well as a list of the order in which the states received seats and which state WOULD have received the next see if one more slot had been available.

CalculateApportionment[total_:435]:= (
	(* Assign one to each state, per the Constitution *)
	stateReps = AssociationThread[stateAbbrs50, ConstantArray[1, 50]];
	(* Keep track of the order in which the remaining (total - 50) seats are assigned by Huntington-Hill *)
	statePriorities = Reverse@Sort@AssociationThread[stateAbbrs50, HuntingtonHill[#, 1]& /@ stateAbbrs50];
	seatOrder = {};
	While[Total@stateReps < total,
		seatOrder = Append[seatOrder, First@Keys@statePriorities];
		{ stateReps, statePriorities } = ApportionOne[stateReps, statePriorities];
	];
	nextSeatWouldBe = First@Keys@statePriorities;
	<| "apportionment" -> stateReps, "order" -> seatOrder, "next" -> nextSeatWouldBe, "size" -> total |>
);

Did it work? Sure did! There difference between the calculated apportionment and the actual one is zero all around

computedApportionment = CalculateApportionment[435];
actualApportionment = AssociationThread[stateAbbrs50, censusData[#]["representatives"]["2010"]& /@ stateAbbrs50];
computedApportionment["apportionment"] - actualApportionment

We’ll want to calculate the ratio of representatives to population. We can’t count DC in the total population since we have no reps! The PeoplePerRepresentative function and VisualizeRepresentation function are available on the Cloud in (HyperlinkURL["visualizeApportionment.nb", https://www.wolframcloud.com/objects/f06fbc86-52d3-4458-afde-d3d5f1bbc27e])

CloudGet["https://www.wolframcloud.com/objects/a14d03ab-e7ec-47b7-90ad-7b1c7a0a2720"]

PeoplePerRepresentation[apportionment, populations]: Calculates each state’s ratio of people to representatives. Called by the visualization function below.

apportionment: Any given 50-state apportionment

populations: The population of each state as an association. Defaults to 2010 decennial Census

VisualizeRepresentation[apportionment, populations, title]: Creates a bar chart of each state’s ratio of people to representatives

apportionment: Any given 50-state apportionment

populations: The population of each state as an association. Defaults to 2010 decennial Census

title: Chart label. Defaults to blank

statePopulations = AssociationThread[
	stateAbbrs50,
	censusData[#]["population_decennial"]["2010"]& /@ stateAbbrs50
];

Total@statePopulations

VisualizeRepresentation[computedApportionment["apportionment"], statePopulations, "Current Apportionment Ratios"]

Let’s see if there’s an ideal size of the House

One of the beauties of the Wolfram language is that it lends itself to running computational experiments. Let’s try a bunch of different House sizes, from 435 members up to three times the current size, and see if there’s an ideal number of representatives. We’ll start by computing all the possible scenarios. This will take just a second since we’re running almost 900 simulations!

houseSizes = Range[435,1305];
allApportionments = AssociationThread[houseSizes, CalculateApportionment /@ houseSizes];

One of the nice things about quantifying “fairness” is that have to define it

There are a few different way we could measure the quality or “fairness” of every scenario. Let’s start by looking at the coefficient of variance (standard deviation / mean) for the people-per-representative values for each scenario.

GetVariance[apportionment_] := (
	ratio = PeoplePerRepresentative[apportionment, statePopulations];
	sd = StandardDeviation[ratio] * 1.0;
	cov = sd / Mean@ratio
)
variances = GetVariance /@ Values@(#["apportionment"]& /@ allApportionments);

How does it look?

ListPlot[Transpose[{houseSizes, variances}], ImageSize->750, AspectRatio->0.3, PlotStyle->{Lighter@Blue, PointSize[0.005]},
		PlotRange->Automatic, PlotRange->All, PlotLabel->Style["Variance by Size of House", FontSize->20], GridLines->{ None, Automatic }]

This is curious!

We’d certainly expect the variance decline as the House grows larger, but there are certain breaking points where one extra seat significantly reduces the disparity, which then climbs up again until the next drop before leveling out. Let’s identify those points before and after these drops to see what’s happening.

indices = Select[Most[Range[Length@variances]], variances[[#]] / variances[[# + 1 ]] > 1.03&] + 1;
optimalSizes = houseSizes[[#]]& /@ indices
optimalValues = variances[[#]]& /@ indices;
optimal = Transpose[{ optimalSizes, optimalValues }];

And the points right before the drop

subOptimalSizes = houseSizes[[#]]& /@ (indices - 1)
subOptimalValues = variances[[#]]& /@ (indices - 1);
subOptimal = Transpose[{ subOptimalSizes, subOptimalValues }];

And plot them to make sure we did that correctly

ListPlot[{Transpose[{houseSizes, variances}], optimal, subOptimal},
	ImageSize -> 900, AspectRatio -> 0.3, PlotRange->Full, Ticks->{True, True},
	PlotStyle->{{Lighter@Blue, Opacity[0.5], PointSize[0.004]},{Green, PointSize[0.009], Opacity[0.75]}, {Red, PointSize[0.009], Opacity[0.75]}},
	PlotLabel->Style["Variance by Size of House, With Key Points", FontSize->20],
	GridLines->{ None, Automatic }
]

Let’s see what the last state to get a representative and the next-in-line state was for each green and red point

summarizeApportionments[apportionments_] := (
	sizes = Flatten[{ Style["Size:", Bold], #["size"]& /@ apportionments }];
	lastSeats = Flatten[{ Style["Last Seat:", Bold], Last@#["order"]& /@ apportionments }];
	lastSeatSizes = Flatten[{ Style["Last Seat #:", Bold], #["apportionment"][Last@#["order"]]& /@ apportionments }];
	nextSeats = Flatten[{ Style["Next-in-line:", Bold], #["next"]& /@ apportionments }];
	nextSeatSizes = Flatten[{ Style["Next-in-line #:", Bold], (#["apportionment"][#["next"]] + 1)& /@ apportionments }];
	Grid[{
		sizes,
		lastSeats,
		lastSeatSizes,
		nextSeats,
		nextSeatSizes
	}, Frame->All, ItemSize->{Automatic, 2}]
)

optimalApportionments = allApportionments[#]& /@ optimalSizes;
summarizeApportionments[optimalApportionments]

subOptimalApportionments = allApportionments[#]& /@ subOptimalSizes;
summarizeApportionments[subOptimalApportionments]

Let’s see what a few of those optimal values look like

populations = AssociationThread[stateAbbrs50, censusData[#]["population_decennial"]["2010"]& /@ stateAbbrs50];

values = { 618, 697, 773, 843, 930, 1207 }; 
ChartValue[value_] := (
	title = ToString@value <> " Seats: Variance of " <> ToString[GetVariance[allApportionments[value]["apportionment"]]];
	Magnify[VisualizeRepresentation[
		CalculateApportionment[value]["apportionment"],
		populations,
		Style[title, FontSize->40]
	], 0.5]
)
Grid[Partition[ChartValue /@ values,2], Frame->All]

Drum roll...

We need 930 members of Congress!

Wait, one last experiment: What if there were 5,000 members?

absurd = CalculateApportionment[5000]["apportionment"];
VisualizeRepresentation[absurd]
GetVariance[absurd]