Repository: econ-ark/HARK (language: python)

Issue:
functional simulation tests to enable external distribution classes
Continuing from this discussion:
https://github.com/econ-ark/HARK/pull/1104#issuecomment-1017620775

I believe that @llorracc is keen on replacing as much of the distribution modeling functionality in HARK with an external library as soon as possible.

A major obstacle to doing this is that the tests for model simulations currently mostly check for near equality with specific numerical values, as opposed to testing for a functional relationship between a model's state value and prior relevant state values.

Fixing this is just a matter of labor time.

This would make it much easier to work on #611 #949 
functional simulation tests to enable external distribution classes
Continuing from this discussion:
https://github.com/econ-ark/HARK/pull/1104#issuecomment-1017620775

I believe that @llorracc is keen on replacing as much of the distribution modeling functionality in HARK with an external library as soon as possible.

A major obstacle to doing this is that the tests for model simulations currently mostly check for near equality with specific numerical values, as opposed to testing for a functional relationship between a model's state value and prior relevant state values.

Fixing this is just a matter of labor time.

This would make it much easier to work on #611 #949 


Discussion / hints:
Though it's not entirely clear from the documentation, I've tested it and this `scipy.stats.rv_discrete` class is more full-featured that our current `DiscreteDistribution` class -- including CDF, entropy, moment, and _expected value_ computation -- with an initialization parameter, `values`, that mirrors the way we are currently constructing DiscreteDistributions.

https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.rv_discrete.html
I've spent a little time in the last couple of days looking into this and learned several things that are relevant for our purposes.  

A key one that illustrates the reason I want to outsource our construction of distributions as much as possible is the recent announcement of [this toolkit](https://arxiv.org/abs/2002.01184) which provides tools for constructing distributions that are adapted to super fast simulation with TensorFlow Markov Chain Monte Carlo tools.  It would be totally unsurprising if someone developed a corresponding tool for PyTorch or other generic tools for doing these kinds of simulations.

Ideally, we could construct an architecture that would make it reasonably straightforward for people to choose `scipy.stats.rv_discrete` or `tfp.mcmc` or the new `StochasticProcess` components of `sympy` when they want to do the Monte Carlo part of their simulations.

Another interesting new development I've just come across is a new toolkit that is laser-focused on describing [bellman](https://github.com/Bellman-devs/bellman) problems, which we should take a very close look at as we think about our future direction.

An example of the kind of test that needs to be changed:

https://github.com/econ-ark/HARK/blob/master/HARK/ConsumptionSaving/tests/test_IndShockConsumerType.py#L92

And example of the kind of test it should be changed into:

https://github.com/econ-ark/HARK/blob/master/HARK/ConsumptionSaving/tests/test_ConsPortfolioFrameModel.py#L153-L156
I have no objection to ADDITIONAL tests like https://github.com/econ-ark/HARK/blob/master/HARK/ConsumptionSaving/tests/test_ConsPortfolioFrameModel.py#L153-L156

And if the reason you think the first test should be eliminated is that it tests the outcome for a specific agent (at index 1), I'm on board.  Instead, what we should be testing is things like whether the average level of assets in simulations matches some raw number like 0.184 (or whatever the correct answer is).

Any time a change in the HARK code causes our simulations to produce significantly different quantitative answers, we need to be alerted to that because we know that our current code produces the right answers and any substantial change in those answers means the HARK change introduced a conceptual error somewhere (even if accounting identities remain satisfied).

Yes, the simulation tests I'm referring to are those that test a specific value for a specific agent.

We also have tests for the solution objects, which should not break with a change in how distributions are sampled. I.e. that the consumption function, for some amount of market resources, is 0.184. So these can go untouched.

I don't recall if we have any tests for the approximate equality of a simulated population mean to a fixed number. I would be wary of such tests, because the sensitivity of the equality approximation would need to be calibrated to the confidence interval of the result somehow to be effective if the only change was to distribution sampling. 
> I don't recall if we have any tests for the approximate equality of a simulated population mean to a fixed number. I would be wary of such tests, because the sensitivity of the equality approximation would need to be calibrated to the confidence interval of the result somehow to be effective if the only change was to distribution sampling.

Any such tests should be low precision ones.  If an aggregate variable changed from 3.14 to 3.15 or 3.13 we would not worry.  But if it changed to 6.8 or something, that would be a clear sign of a bug.

As of the meeting today, we decided:
 - For numerical tests where precision is an issue, change the precision (not in the scope of this ticket)
 - For numerical tests of simulation outputs, which are blocking some changes moving forward, we can remove those current tests (!)
   - This reduces test coverage. Functional simulation tests as discussed above are still a good idea so that there is test coverage of the transition equations. 
I'm taking responsibility for this issue, but am a little stuck.

One thing we could do is remove any tests for particular values that come out of a simulation's results.

In many cases, this would leave simulations as entirely untested except to show that they executed without error.

I wonder if there's anything else that comes to mind for what would be a good simulation test.
(The tests of a particular transition function are one idea, though writing these will be a cumbersome process.)
Though it's not entirely clear from the documentation, I've tested it and this `scipy.stats.rv_discrete` class is more full-featured that our current `DiscreteDistribution` class -- including CDF, entropy, moment, and _expected value_ computation -- with an initialization parameter, `values`, that mirrors the way we are currently constructing DiscreteDistributions.

https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.rv_discrete.html
I've spent a little time in the last couple of days looking into this and learned several things that are relevant for our purposes.  

A key one that illustrates the reason I want to outsource our construction of distributions as much as possible is the recent announcement of [this toolkit](https://arxiv.org/abs/2002.01184) which provides tools for constructing distributions that are adapted to super fast simulation with TensorFlow Markov Chain Monte Carlo tools.  It would be totally unsurprising if someone developed a corresponding tool for PyTorch or other generic tools for doing these kinds of simulations.

Ideally, we could construct an architecture that would make it reasonably straightforward for people to choose `scipy.stats.rv_discrete` or `tfp.mcmc` or the new `StochasticProcess` components of `sympy` when they want to do the Monte Carlo part of their simulations.

Another interesting new development I've just come across is a new toolkit that is laser-focused on describing [bellman](https://github.com/Bellman-devs/bellman) problems, which we should take a very close look at as we think about our future direction.

An example of the kind of test that needs to be changed:

https://github.com/econ-ark/HARK/blob/master/HARK/ConsumptionSaving/tests/test_IndShockConsumerType.py#L92

And example of the kind of test it should be changed into:

https://github.com/econ-ark/HARK/blob/master/HARK/ConsumptionSaving/tests/test_ConsPortfolioFrameModel.py#L153-L156
I have no objection to ADDITIONAL tests like https://github.com/econ-ark/HARK/blob/master/HARK/ConsumptionSaving/tests/test_ConsPortfolioFrameModel.py#L153-L156

And if the reason you think the first test should be eliminated is that it tests the outcome for a specific agent (at index 1), I'm on board.  Instead, what we should be testing is things like whether the average level of assets in simulations matches some raw number like 0.184 (or whatever the correct answer is).

Any time a change in the HARK code causes our simulations to produce significantly different quantitative answers, we need to be alerted to that because we know that our current code produces the right answers and any substantial change in those answers means the HARK change introduced a conceptual error somewhere (even if accounting identities remain satisfied).

Yes, the simulation tests I'm referring to are those that test a specific value for a specific agent.

We also have tests for the solution objects, which should not break with a change in how distributions are sampled. I.e. that the consumption function, for some amount of market resources, is 0.184. So these can go untouched.

I don't recall if we have any tests for the approximate equality of a simulated population mean to a fixed number. I would be wary of such tests, because the sensitivity of the equality approximation would need to be calibrated to the confidence interval of the result somehow to be effective if the only change was to distribution sampling. 
> I don't recall if we have any tests for the approximate equality of a simulated population mean to a fixed number. I would be wary of such tests, because the sensitivity of the equality approximation would need to be calibrated to the confidence interval of the result somehow to be effective if the only change was to distribution sampling.

Any such tests should be low precision ones.  If an aggregate variable changed from 3.14 to 3.15 or 3.13 we would not worry.  But if it changed to 6.8 or something, that would be a clear sign of a bug.

As of the meeting today, we decided:
 - For numerical tests where precision is an issue, change the precision (not in the scope of this ticket)
 - For numerical tests of simulation outputs, which are blocking some changes moving forward, we can remove those current tests (!)
   - This reduces test coverage. Functional simulation tests as discussed above are still a good idea so that there is test coverage of the transition equations. 
I'm taking responsibility for this issue, but am a little stuck.

One thing we could do is remove any tests for particular values that come out of a simulation's results.

In many cases, this would leave simulations as entirely untested except to show that they executed without error.

I wonder if there's anything else that comes to mind for what would be a good simulation test.
(The tests of a particular transition function are one idea, though writing these will be a cumbersome process.)

Fix patch:
diff --git a/Documentation/CHANGELOG.md b/Documentation/CHANGELOG.md
index b0c847f26..992a38d68 100644
--- a/Documentation/CHANGELOG.md
+++ b/Documentation/CHANGELOG.md
@@ -40,6 +40,7 @@ Release Date: TBD
 * Removes the extra-dimension that was returned by `calc_expectations` in some instances. [#1149](https://github.com/econ-ark/HARK/pull/1149)
 * Adds `HARK.distribution.expected` alias for `DiscreteDistribution.expected`. [#1156](https://github.com/econ-ark/HARK/pull/1156)
 * Renames attributes in `DiscreteDistribution`: `X` to `atoms` and `pmf` to `pmv`. [#1164](https://github.com/econ-ark/HARK/pull/1164), [#1051](https://github.com/econ-ark/HARK/pull/1151), [#1159](https://github.com/econ-ark/HARK/pull/1159). 
+* Remove or replace automated tests that depend on brittle simulation results. [#1148](https://github.com/econ-ark/HARK/pull/1148)
 * Updates asset grid constructor from `ConsIndShockModel.py` to allow for linearly-spaced grids when `aXtraNestFac == -1`. [#1172](https://github.com/econ-ark/HARK/pull/1172)
 
 ### 0.12.0


Test patch:
diff --git a/HARK/ConsumptionSaving/tests/test_ConsAggShockModel.py b/HARK/ConsumptionSaving/tests/test_ConsAggShockModel.py
index 2eb2be79a..fedaa51bb 100644
--- a/HARK/ConsumptionSaving/tests/test_ConsAggShockModel.py
+++ b/HARK/ConsumptionSaving/tests/test_ConsAggShockModel.py
@@ -54,7 +54,8 @@ def test_macro(self):
             self.economy.AFunc.slope, 1.116259456228145
         )
 
-        self.assertAlmostEqual(self.economy.history["MaggNow"][10], 7.456324335623432)
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(self.economy.history["MaggNow"][10], 7.456324335623432)
 
 
 class testAggShockMarkovConsumerType(unittest.TestCase):
@@ -90,7 +91,8 @@ def test_economy(self):
             self.economy.AFunc[0].slope, 1.0904698841958917
         )
 
-        self.assertAlmostEqual(self.economy.history["AaggNow"][5], 9.467758924955897)
+        # simulation test -- seed/generator specific
+        # self.assertAlmostEqual(self.economy.history["AaggNow"][5], 9.467758924955897)
 
 
 class KrusellSmithTestCase(unittest.TestCase):
@@ -150,10 +152,12 @@ def test_methods(self):
         )
 
         self.economy.make_Mrkv_history()  # Make a simulated history of aggregate shocks
-        self.assertAlmostEqual(
-            np.mean(self.economy.MrkvNow_hist),
-            0.4818181818181818
-        )
+
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(
+        #    np.mean(self.economy.MrkvNow_hist),
+        #    0.4818181818181818
+        #)
 
         # object attributes that are conditions
         # for precompute_arrays
@@ -185,15 +189,16 @@ def test_methods(self):
 
         # testing update_solution_terminal()
         self.assertEqual(
-            self.agent.solution_terminal.cFunc[0](10,self.economy.MSS),
+            self.agent.solution_terminal.cFunc[0](10, self.economy.MSS),
             10
         )
 
         self.assertAlmostEqual(
             self.economy.agents[0].solution[0].cFunc[0](
-                10,self.economy.MSS
+                10, self.economy.MSS
             ).tolist(),
-            1.2386775112633517
+            1.2386774,
+            places = 4
         )
 
         self.assertAlmostEqual(
@@ -205,29 +210,34 @@ def test_methods(self):
 
         emp_totals = np.sum(self.agent.history['EmpNow'], axis = 0)
 
-        self.assertEqual(emp_totals[0], 1011)
-        self.assertEqual(emp_totals[2], 1009)
-        self.assertEqual(emp_totals[9], 1042)
-
-        self.assertAlmostEqual(
-            self.economy.history['Aprev'][0],
-            11.83133152125537
-        )
-
-        self.assertAlmostEqual(
-            self.economy.history['Aprev'][1],
-            11.260758949042032
-        )
-
-        self.assertAlmostEqual(
-            self.economy.history['Aprev'][2],
-            10.72309150585389
-        )
-
-        self.assertAlmostEqual(
-            self.economy.history['Mnow'][10],
-            self.economy.history['Mnow'][10]
-        )
+        # simulation test -- seed/generator specific
+        # self.assertEqual(emp_totals[0], 1011)
+        # self.assertEqual(emp_totals[2], 1009)
+        # self.assertEqual(emp_totals[9], 1042)
+
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(
+        #    self.economy.history['Aprev'][0],
+        #    11.83133152125537
+        #)
+
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(
+        #    self.economy.history['Aprev'][1],
+        #    11.260758949042032
+        #)
+
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(
+        #    self.economy.history['Aprev'][2],
+        #    10.72309150585389
+        #)
+
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(
+        #    self.economy.history['Mnow'][10],
+        #    self.economy.history['Mnow'][10]
+        #)
 
         new_dynamics = self.economy.update_dynamics()
 
@@ -261,8 +271,11 @@ def test_economy(self):
         self.economy.AFunc = self.economy.dynamics.AFunc
         self.assertAlmostEqual(self.economy.AFunc[0].slope, 1.0014463644834282)
 
-        self.assertAlmostEqual(self.economy.history["Aprev"][4], 11.009107526443584)
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(self.economy.history["Aprev"][4], 11.009107526443584)
 
-        self.assertAlmostEqual(self.economy.history['Mrkv'][40], 1)
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(self.economy.history['Mrkv'][40], 1)
 
-        self.assertAlmostEqual(self.economy.history["Urate"][12], 0.040000000000000036)
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(self.economy.history["Urate"][12], 0.040000000000000036)
diff --git a/HARK/ConsumptionSaving/tests/test_ConsGenIncProcessModel.py b/HARK/ConsumptionSaving/tests/test_ConsGenIncProcessModel.py
index 4602a057b..679d72fa2 100644
--- a/HARK/ConsumptionSaving/tests/test_ConsGenIncProcessModel.py
+++ b/HARK/ConsumptionSaving/tests/test_ConsGenIncProcessModel.py
@@ -100,4 +100,5 @@ def test_simulation(self):
         self.agent.initialize_sim()
         self.agent.simulate()
 
-        self.assertAlmostEqual(np.mean(self.agent.history["mLvl"]), 1.2043946738813716)
+        # simulation test -- seed/generator specific
+        # self.assertAlmostEqual(np.mean(self.agent.history["mLvl"]), 1.2043946738813716)
diff --git a/HARK/ConsumptionSaving/tests/test_ConsPortfolioFrameModel.py b/HARK/ConsumptionSaving/tests/test_ConsPortfolioFrameModel.py
index 37559d5a1..ef51f72e2 100644
--- a/HARK/ConsumptionSaving/tests/test_ConsPortfolioFrameModel.py
+++ b/HARK/ConsumptionSaving/tests/test_ConsPortfolioFrameModel.py
@@ -36,15 +36,19 @@ def test_simOnePeriod(self):
         self.pcct.track_vars += ["aNrm"]
         self.pcct.initialize_sim()
 
-        self.assertFalse(np.any(self.pcct.shocks["Adjust"]))
+        # simulation test -- seed/generator specific
+        # self.assertFalse(np.any(self.pcct.shocks["Adjust"]))
 
         self.pcct.sim_one_period()
 
-        self.assertAlmostEqual(self.pcct.shocks["PermShk"][0], 0.9692322)
+        # simulation test -- seed/generator specific
+        # self.assertAlmostEqual(self.pcct.shocks["PermShk"][0], 0.9692322)
 
-        self.assertAlmostEqual(self.pcct.shocks["TranShk"][0], 1.03172631)
+        # simulation test -- seed/generator specific
+        # self.assertAlmostEqual(self.pcct.shocks["TranShk"][0], 1.03172631)
 
-        self.assertAlmostEqual(self.pcct.shocks["Risky"][0], 0.96358739)
+        # simulation test -- seed/generator specific
+        # self.assertAlmostEqual(self.pcct.shocks["Risky"][0], 0.96358739)
 
         self.assertAlmostEqual(
             self.pcct.state_now["pLvl"][0],
@@ -61,11 +65,13 @@ def test_simOnePeriod(self):
             + self.pcct.shocks["TranShk"][0],
         )
 
-        self.assertAlmostEqual(
-            # todo: more flexible test
-            self.pcct.controls["Share"][0],
-            0.90256316,
-        )
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(
+        #    # todo: more flexible test
+        #    self.pcct.controls["Share"][0],
+        #    0.90256316,
+        #)
+        
         self.assertAlmostEqual(
             self.pcct.controls["cNrm"][0],
             self.pcct.solution[0].cFuncAdj(self.pcct.state_now["mNrm"][0]),
diff --git a/HARK/ConsumptionSaving/tests/test_ConsPortfolioModel.py b/HARK/ConsumptionSaving/tests/test_ConsPortfolioModel.py
index 089d98452..c09dbbb20 100644
--- a/HARK/ConsumptionSaving/tests/test_ConsPortfolioModel.py
+++ b/HARK/ConsumptionSaving/tests/test_ConsPortfolioModel.py
@@ -39,37 +39,47 @@ def test_solution(self):
     def test_sim_one_period(self):
         self.pcct.T_sim = 30
         self.pcct.AgentCount = 10
-        self.pcct.track_vars += ["aNrm"]
+        self.pcct.track_vars += ["aNrm", "mNrm", "bNrm", "TranShk", "cNrm"]
         self.pcct.initialize_sim()
 
         self.assertFalse(np.any(self.pcct.shocks["Adjust"]))
 
-        self.assertAlmostEqual(self.pcct.state_now["pLvl"][0], 1.0)
+        # simulation test -- seed/generator specific
+        # self.assertAlmostEqual(self.pcct.state_now["pLvl"][0], 1.0)
 
-        self.assertAlmostEqual(self.pcct.state_now["aNrm"][0], 7.257027956)
+        # simulation test -- seed/generator specific
+        # self.assertAlmostEqual(self.pcct.state_now["aNrm"][0], 7.257027956)
 
-        self.assertAlmostEqual(self.pcct.Rfree[0], 1.03)
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(self.pcct.Rfree[0], 1.03)
 
-        self.assertAlmostEqual(self.pcct.state_now["PlvlAgg"], 1.0)
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(self.pcct.state_now["PlvlAgg"], 1.0)
 
         self.pcct.sim_one_period()
 
-        self.assertAlmostEqual(self.pcct.state_now["pLvl"][0], 0.858934461)
-
-        self.assertAlmostEqual(self.pcct.state_now["aNrm"][0], 8.02359093)
-
-        self.assertAlmostEqual(self.pcct.shocks["PermShk"][0], 0.85893446)
-
-        self.assertAlmostEqual(self.pcct.shocks["TranShk"][0], 1.0)
+        self.assertAlmostEqual(
+            self.pcct.state_now["mNrm"][0],
+            self.pcct.state_now["bNrm"][0] + self.pcct.shocks["TranShk"][0],
+        )
 
-        self.assertAlmostEqual(self.pcct.shocks["Risky"][0], 0.92405816)
+        self.assertAlmostEqual(
+            self.pcct.controls["cNrm"][0],
+            self.pcct.solution[0].cFuncAdj(self.pcct.state_now["mNrm"][0]),
+        )
 
-        self.assertTrue(np.any(self.pcct.shocks["Adjust"][0]))
+        self.assertAlmostEqual(
+            self.pcct.controls["Share"][0],
+            self.pcct.solution[0].ShareFuncAdj(self.pcct.state_now["mNrm"][0]),
+        )
 
-        self.assertAlmostEqual(self.pcct.state_now["mNrm"][0], 9.70233892039)
+        self.assertAlmostEqual(
+            self.pcct.state_now["aNrm"][0],
+            self.pcct.state_now["mNrm"][0] - self.pcct.controls["cNrm"][0],
+        )
 
-        self.assertAlmostEqual(self.pcct.controls["Share"][0], 0.8627164488246847)
-        self.assertAlmostEqual(self.pcct.controls["cNrm"][0], 1.67874799)
+        # a drawn shock ; may not be robust to RNG/disitrubition implementations
+        #self.assertAlmostEqual(self.pcct.shocks["Adjust"][0], 1.0)
 
 
 class SimulatePortfolioConsumerTypeTestCase(PortfolioConsumerTypeTestCase):
@@ -86,42 +96,60 @@ def test_simulation(self):
             "Adjust",
             "PermShk",
             "bNrm",
+            "TranShk"
         ]
         self.pcct.initialize_sim()
 
         self.pcct.simulate()
 
-        self.assertAlmostEqual(self.pcct.history["mNrm"][0][0], 9.70233892)
+        self.assertAlmostEqual(
+            self.pcct.history["mNrm"][0][0],
+            self.pcct.history["bNrm"][0][0] + self.pcct.history["TranShk"][0][0],
+        )
 
-        self.assertAlmostEqual(self.pcct.history["cNrm"][0][0], 1.6787479894848298)
+        self.assertAlmostEqual(
+            self.pcct.history["cNrm"][0][0],
+            self.pcct.solution[0].cFuncAdj(self.pcct.history["mNrm"][0][0]),
+        )
 
-        self.assertAlmostEqual(self.pcct.history["Share"][0][0], 0.8627164488246847)
+        self.assertAlmostEqual(
+            self.pcct.history["Share"][0][0],
+            self.pcct.solution[0].ShareFuncAdj(self.pcct.history["mNrm"][0][0]),
+        )
 
-        self.assertAlmostEqual(self.pcct.history["aNrm"][0][0], 8.023590930905383)
+        self.assertAlmostEqual(
+            self.pcct.history["aNrm"][0][0],
+            self.pcct.history["mNrm"][0][0] - self.pcct.history["cNrm"][0][0],
+        )
 
         self.assertAlmostEqual(self.pcct.history["Adjust"][0][0], 1.0)
 
         # the next period
-        self.assertAlmostEqual(self.pcct.history["Risky"][1][0], 0.8950304697526602)
-
-        self.assertAlmostEqual(self.pcct.history["Rport"][1][0], 0.9135595661654792)
-
-        self.assertAlmostEqual(self.pcct.history["Adjust"][1][0], 1.0)
-
-        self.assertAlmostEqual(self.pcct.history["PermShk"][1][0], 1.0050166461586711)
 
-        self.assertAlmostEqual(self.pcct.history["bNrm"][1][0], 7.293439643953855)
-
-        self.assertAlmostEqual(self.pcct.history["mNrm"][1][0], 8.287859049575047)
-
-        self.assertAlmostEqual(self.pcct.history["cNrm"][1][0], 1.5773607434989751)
+        self.assertAlmostEqual(
+            self.pcct.history["mNrm"][1][0],
+            self.pcct.history["bNrm"][1][0] + self.pcct.history["TranShk"][1][0],
+        )
 
-        self.assertAlmostEqual(self.pcct.history["Share"][1][0], 0.9337608822146805)
+        self.assertAlmostEqual(
+            self.pcct.history["cNrm"][1][0],
+            self.pcct.solution[0].cFuncAdj(self.pcct.history["mNrm"][1][0]),
+        )
 
-        self.assertAlmostEqual(self.pcct.history["aNrm"][1][0], 6.710498306076072)
+        self.assertAlmostEqual(
+            self.pcct.history["Share"][1][0],
+            self.pcct.solution[0].ShareFuncAdj(self.pcct.history["mNrm"][1][0]),
+        )
 
-        self.assertAlmostEqual(self.pcct.history["aNrm"][15][0], 5.304746367434934)
+        self.assertAlmostEqual(
+            self.pcct.history["aNrm"][1][0],
+            self.pcct.history["mNrm"][1][0] - self.pcct.history["cNrm"][1][0],
+        )
 
+        self.assertAlmostEqual(
+            self.pcct.history["aNrm"][15][0],
+            self.pcct.history["mNrm"][15][0] - self.pcct.history["cNrm"][15][0],
+        )
 
 class testPortfolioConsumerTypeSticky(unittest.TestCase):
     def test_sticky(self):
diff --git a/HARK/ConsumptionSaving/tests/test_ConsPrefShockModel.py b/HARK/ConsumptionSaving/tests/test_ConsPrefShockModel.py
index eeb44ca6e..c676c6d91 100644
--- a/HARK/ConsumptionSaving/tests/test_ConsPrefShockModel.py
+++ b/HARK/ConsumptionSaving/tests/test_ConsPrefShockModel.py
@@ -36,14 +36,16 @@ def test_simulation(self):
         self.agent.initialize_sim()
         self.agent.simulate()
 
-        self.assertAlmostEqual(self.agent.history["cNrm"][0][5], 0.7366020536567589)
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(self.agent.history["cNrm"][0][5], 0.7366020536567589)
 
         self.assertEqual(
             self.agent.shock_history["PrefShk"][0][5],
             self.agent.history["PrefShk"][0][5],
         )
 
-        self.assertEqual(self.agent.history["PrefShk"][0][5], 0.4909415933881665)
+        # simulation test -- seed/generator specific
+        # self.assertEqual(self.agent.history["PrefShk"][0][5], 0.4909415933881665)
 
 
 class testKinkyPrefConsumerType(unittest.TestCase):
@@ -77,4 +79,5 @@ def test_simulation(self):
         self.agent.initialize_sim()
         self.agent.simulate()
 
-        self.assertAlmostEqual(self.agent.history["cNrm"][0][5], 0.7717096928111515)
+        # simulation test -- seed/generator specific
+        # self.assertAlmostEqual(self.agent.history["cNrm"][0][5], 0.7717096928111515)
diff --git a/HARK/ConsumptionSaving/tests/test_IndShockConsumerType.py b/HARK/ConsumptionSaving/tests/test_IndShockConsumerType.py
index 1e6eec257..fe78daf6d 100644
--- a/HARK/ConsumptionSaving/tests/test_IndShockConsumerType.py
+++ b/HARK/ConsumptionSaving/tests/test_IndShockConsumerType.py
@@ -25,9 +25,10 @@ def test_get_shocks(self):
 
         self.agent.get_shocks()
 
-        self.assertEqual(self.agent.shocks["PermShk"][0], 1.0427376294215103)
-        self.assertAlmostEqual(self.agent.shocks["PermShk"][1], 0.9278094171517413)
-        self.assertAlmostEqual(self.agent.shocks["TranShk"][0], 0.881761797501595)
+        # simulation test -- seed/generator specific
+        #self.assertEqual(self.agent.shocks["PermShk"][0], 1.0427376294215103)
+        #self.assertAlmostEqual(self.agent.shocks["PermShk"][1], 0.9278094171517413)
+        #self.assertAlmostEqual(self.agent.shocks["TranShk"][0], 0.881761797501595)
 
     def test_ConsIndShockSolverBasic(self):
         LifecycleExample = IndShockConsumerType(**init_lifecycle)
@@ -93,9 +94,11 @@ def test_simulated_values(self):
         self.agent.initialize_sim()
         self.agent.simulate()
 
-        self.assertAlmostEqual(self.agent.MPCnow[1], 0.5711503906043797)
-
-        self.assertAlmostEqual(self.agent.state_now["aLvl"][1], 0.18438326264597635)
+        # MPCnow depends on assets, which are stochastic
+        #self.assertAlmostEqual(self.agent.MPCnow[1], 0.5711503906043797)
+        
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(self.agent.state_now["aLvl"][1], 0.18438326264597635)
 
 
 class testBufferStock(unittest.TestCase):
@@ -248,9 +251,10 @@ def test_infinite_horizon(self):
         IndShockExample.initialize_sim()
         IndShockExample.simulate()
 
-        self.assertAlmostEqual(
-            IndShockExample.history["mNrm"][0][0], 1.0170176090252379
-        )
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(
+        #    IndShockExample.history["mNrm"][0][0], 1.0170176090252379
+        #)
 
 
 LifecycleDict = {  # Click arrow to expand this fairly large parameter dictionary
@@ -698,9 +702,12 @@ def test_Harmenberg_mtd(self):
         c_std1 = np.std(Consumption_list)
         c_std_ratio = c_std2 / c_std1
 
-        self.assertAlmostEqual(c_std2, 0.03768819564871894)
-        self.assertAlmostEqual(c_std1, 0.004411745897568616)
-        self.assertAlmostEqual(c_std_ratio, 8.542694099741672)
+        # simulation tests -- seed/generator specific
+        # But these are based on aggregate population statistics.
+        # WARNING: May fail stochastically, or based on specific RNG types.
+        #self.assertAlmostEqual(c_std2, 0.0376882, places = 2)
+        #self.assertAlmostEqual(c_std1, 0.0044117, places = 2)
+        #self.assertAlmostEqual(c_std_ratio, 8.5426941, places = 2)
 
 
 # %% Shock pre-computing tests
diff --git a/HARK/ConsumptionSaving/tests/test_IndShockConsumerTypeFast.py b/HARK/ConsumptionSaving/tests/test_IndShockConsumerTypeFast.py
index ae3da7485..481795e57 100644
--- a/HARK/ConsumptionSaving/tests/test_IndShockConsumerTypeFast.py
+++ b/HARK/ConsumptionSaving/tests/test_IndShockConsumerTypeFast.py
@@ -27,9 +27,10 @@ def test_get_shocks(self):
 
         self.agent.get_shocks()
 
-        self.assertEqual(self.agent.shocks['PermShk'][0], 1.0427376294215103)
-        self.assertAlmostEqual(self.agent.shocks['PermShk'][1], 0.9278094171517413)
-        self.assertAlmostEqual(self.agent.shocks['TranShk'][0], 0.881761797501595)
+        # simulation test -- seed/generator specific
+        #self.assertEqual(self.agent.shocks['PermShk'][0], 1.0427376294215103)
+        #self.assertAlmostEqual(self.agent.shocks['PermShk'][1], 0.9278094171517413)
+        #self.assertAlmostEqual(self.agent.shocks['TranShk'][0], 0.881761797501595)
 
     def test_ConsIndShockSolverBasic(self):
         LifecycleExample = IndShockConsumerTypeFast(**init_lifecycle)
@@ -56,9 +57,12 @@ def test_ConsIndShockSolverBasic(self):
     def test_simulated_values(self):
         self.agent.initialize_sim()
         self.agent.simulate()
-        self.assertAlmostEqual(self.agent.MPCnow[1], 0.5711503906043797)
+        
+        #MPCnow is stochastic
+        #self.assertAlmostEqual(self.agent.MPCnow[1], 0.5711503906043797)
 
-        self.assertAlmostEqual(self.agent.state_now['aLvl'][1], 0.18438326264597635)
+        # simulation test -- seed/generator specific
+        # self.assertAlmostEqual(self.agent.state_now['aLvl'][1], 0.18438326264597635)
 
 
 class testBufferStock(unittest.TestCase):
@@ -169,9 +173,10 @@ def test_infinite_horizon(self):
         IndShockExample.initialize_sim()
         IndShockExample.simulate()
 
-        self.assertAlmostEqual(
-            IndShockExample.history["mNrm"][0][0], 1.0170176090252379
-        )
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(
+        #    IndShockExample.history["mNrm"][0][0], 1.0170176090252379
+        #)
 
 
 class testIndShockConsumerTypeFastLifecycle(unittest.TestCase):
diff --git a/HARK/ConsumptionSaving/tests/test_PerfForesightConsumerType.py b/HARK/ConsumptionSaving/tests/test_PerfForesightConsumerType.py
index 1ae385ac3..c5f4d1424 100644
--- a/HARK/ConsumptionSaving/tests/test_PerfForesightConsumerType.py
+++ b/HARK/ConsumptionSaving/tests/test_PerfForesightConsumerType.py
@@ -65,19 +65,20 @@ def test_simulation(self):
         )  # This implicitly uses the assign_parameters method of AgentType
 
         # Create PFexample object
-        self.agent_infinite.track_vars = ["mNrm"]
+        self.agent_infinite.track_vars = ["bNrm", "mNrm", "TranShk"]
         self.agent_infinite.initialize_sim()
         self.agent_infinite.simulate()
 
         self.assertAlmostEqual(
             np.mean(self.agent_infinite.history["mNrm"], axis=1)[40],
-            -23.008063500363942,
+            np.mean(self.agent_infinite.history["bNrm"], axis=1)[40] + np.mean(self.agent_infinite.history["TranShk"], axis=1)[40]
         )
 
-        self.assertAlmostEqual(
-            np.mean(self.agent_infinite.history["mNrm"], axis=1)[100],
-            -27.164608851546927,
-        )
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(
+        #    np.mean(self.agent_infinite.history["mNrm"], axis=1)[100],
+        #    -27.164608851546927,
+        #)
 
         ## Try now with the manipulation at time step 80
 
@@ -91,15 +92,17 @@ def test_simulation(self):
         )
         self.agent_infinite.simulate(40)
 
-        self.assertAlmostEqual(
-            np.mean(self.agent_infinite.history["mNrm"], axis=1)[40],
-            -23.008063500363942,
-        )
-
-        self.assertAlmostEqual(
-            np.mean(self.agent_infinite.history["mNrm"], axis=1)[100],
-            -29.140261331951606,
-        )
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(
+        #    np.mean(self.agent_infinite.history["mNrm"], axis=1)[40],
+        #    -23.008063500363942,
+        #)
+
+        # simulation test -- seed/generator specific
+        #self.assertAlmostEqual(
+        #    np.mean(self.agent_infinite.history["mNrm"], axis=1)[100],
+        #    -29.140261331951606,
+        #)
         
     def test_stable_points(self):
         
diff --git a/HARK/ConsumptionSaving/tests/test_TractableBufferStockModel.py b/HARK/ConsumptionSaving/tests/test_TractableBufferStockModel.py
index 0c0cfcc52..6d8b45a4d 100644
--- a/HARK/ConsumptionSaving/tests/test_TractableBufferStockModel.py
+++ b/HARK/ConsumptionSaving/tests/test_TractableBufferStockModel.py
@@ -16,8 +16,8 @@ def test_simulation(self):
         self.tct.aLvlInitMean = 0.0  # Mean of log initial assets for new agents
         self.tct.aLvlInitStd = 1.0 # stdev of log initial assets for new agents
         self.tct.T_cycle = 1
-        self.tct.track_vars += ["mLvl"]
+        self.tct.track_vars += ["mLvl", "cLvlNow", "aLvl"]
         self.tct.initialize_sim()
         self.tct.simulate()
 
-        self.assertAlmostEqual(self.tct.history["mLvl"][15][0], 5.820630251772332)
+        self.assertAlmostEqual(self.tct.history["mLvl"][15][0] - self.tct.history["cLvlNow"][15][0], self.tct.history["aLvl"][15][0] )
diff --git a/HARK/tests/test_distribution.py b/HARK/tests/test_distribution.py
index 7d56ea1eb..7887c5134 100644
--- a/HARK/tests/test_distribution.py
+++ b/HARK/tests/test_distribution.py
@@ -285,54 +285,56 @@ class DistributionClassTests(unittest.TestCase):
     """
 
     def test_drawMeanOneLognormal(self):
-        self.assertEqual(MeanOneLogNormal().draw(1)[0], 3.5397367004222002)
+        MeanOneLogNormal().draw(1)[0]
 
     def test_Lognormal(self):
         dist = Lognormal()
 
-        self.assertEqual(dist.draw(1)[0], 5.836039190663969)
+        dist.draw(1)[0]
 
         dist.draw(100)
         dist.reset()
 
-        self.assertEqual(dist.draw(1)[0], 5.836039190663969)
+        dist.draw(1)[0]
 
     def test_Normal(self):
         dist = Normal()
 
-        self.assertEqual(dist.draw(1)[0], 1.764052345967664)
+        dist.draw(1)[0]
 
         dist.draw(100)
         dist.reset()
 
-        self.assertEqual(dist.draw(1)[0], 1.764052345967664)
+        dist.draw(1)[0]
 
     def test_MVNormal(self):
+
+        ## Are these tests generator/backend specific?
         dist = MVNormal()
 
-        self.assertTrue(
-            np.allclose(dist.draw(1)[0], np.array([2.76405235, 1.40015721]))
-        )
+        #self.assertTrue(
+        #    np.allclose(dist.draw(1)[0], np.array([2.76405235, 1.40015721]))
+        #)
 
         dist.draw(100)
         dist.reset()
 
-        self.assertTrue(
-            np.allclose(dist.draw(1)[0], np.array([2.76405235, 1.40015721]))
-        )
+        #self.assertTrue(
+        #    np.allclose(dist.draw(1)[0], np.array([2.76405235, 1.40015721]))
+        #)
 
     def test_Weibull(self):
-        self.assertAlmostEqual(Weibull().draw(1)[0], 0.79587450816311)
+        Weibull().draw(1)[0]
 
     def test_Uniform(self):
         uni = Uniform()
 
-        self.assertEqual(Uniform().draw(1)[0], 0.5488135039273248)
+        Uniform().draw(1)[0]
 
         self.assertEqual(calc_expectation(uni.approx(10)), 0.5)
 
     def test_Bernoulli(self):
-        self.assertEqual(Bernoulli().draw(1)[0], False)
+        Bernoulli().draw(1)[0]
 
 
 class IndexDistributionClassTests(unittest.TestCase):