diff --git a/ogcore/aggregates.py b/ogcore/aggregates.py index 63b2a94a6..d7d3393ce 100644 --- a/ogcore/aggregates.py +++ b/ogcore/aggregates.py @@ -133,8 +133,8 @@ def get_B(b, p, method, preTP): if preTP: part1 = b * np.transpose(p.omega_S_preTP * p.lambdas) omega_extended = np.append(p.omega_S_preTP[1:], [0.0]) - imm_extended = np.append(p.imm_rates[0, 1:], [0.0]) - pop_growth_rate = p.g_n[0] + imm_extended = np.append(p.imm_rates_preTP[1:], [0.0]) + pop_growth_rate = p.g_n_preTP else: part1 = b * np.transpose(p.omega_SS * p.lambdas) omega_extended = np.append(p.omega_SS[1:], [0.0]) @@ -189,8 +189,8 @@ def get_BQ(r, b_splus1, j, p, method, preTP): if method == "SS": if preTP: omega = p.omega_S_preTP - pop_growth_rate = p.g_n[0] - rho = p.rho[0, :] + pop_growth_rate = p.g_n_preTP + rho = p.rho_preTP else: omega = p.omega_SS pop_growth_rate = p.g_n_ss @@ -206,7 +206,7 @@ def get_BQ(r, b_splus1, j, p, method, preTP): p.omega_S_preTP.reshape(1, p.S), p.omega[: p.T - 1, :], axis=0 ) rho = np.append( - p.rho[0, :].reshape(1, p.S), p.rho[: p.T - 1, :], axis=0 + p.rho_preTP.reshape(1, p.S), p.rho[: p.T - 1, :], axis=0 ) if j is not None: diff --git a/ogcore/default_parameters.json b/ogcore/default_parameters.json index 2922d15f3..3bd62b0fe 100644 --- a/ogcore/default_parameters.json +++ b/ogcore/default_parameters.json @@ -4531,8 +4531,7 @@ "value": [ { "value": [0.017479913702712725, 0.017252673721361574, 0.017313468985451712, 0.01741734467720022, 0.017610836147653894, 0.018037326972746067, 0.018428614504312973, 0.018731212768297128, 0.019223877430412158, 0.019533665622202632, 0.019434917008899554, 0.018636460788059497, 0.018161740350649384, 0.017819885884360336, 0.017856106761386912, 0.017939064896242172, 0.0173584158450137, 0.01759326731381578, 0.017556092457882918, 0.017348831455567988, 0.01763610875892841, 0.016521134511611055, 0.01616493583053726, 0.015964729480990054, 0.015477015195103342, 0.01591654299984056, 0.01539650128360123, 0.01561223197188071, 0.016292741432876165, 0.017193018202081137, 0.01736681745887612, 0.016428597852015135, 0.01605376707200027, 0.016036840400303197, 0.016306172531601798, 0.017250621809611313, 0.017536449193224023, 0.017487907335860874, 0.017416472477538542, 0.017590642752769312, 0.01769271828193883, 0.017056672755048446, 0.016903326974848074, 0.016655018441783655, 0.01601349080436947, 0.015817814264164037, 0.015095305117520927, 0.014438237313875455, 0.013801453267642878, 0.013298311958716617, 0.012889796775514745, 0.012405992316545383, 0.012171305937986552, 0.012537090181191982, 0.009140572473656714, 0.008895104621663985, 0.00854997250414294, 0.008604276207804186, 0.007359240011714394, 0.0065857742259522725, 0.00608824742682594, 0.005584969530483296, 0.005194883294310685, 0.00463722646525393, 0.004268302269151332, 0.003921275687879491, 0.003360617443596329, 0.00307989603033834, 0.0027689270722019797, 0.002464564334019259, 0.002189389475001142, 0.001817465282607566, 0.0015516186002795162, 0.001280287314428285, 0.0010211214995148833, 0.0008163607706652746, 0.0006238227539144882, 0.0004616610853196622, 0.0003378106654761438, 0.00023711271717107317] - } - ], + } ], "validators": { "range": { "min": 0.0, @@ -4577,6 +4576,24 @@ } } }, + "g_n_preTP": { + "title": "Population growth rate from year before model start to start year", + "description": "Population growth rate from year before model start to start year.", + "section_1": "Demographic Parameters", + "notes": "", + "type": "float", + "value": [ + { + "value": 0.0074789353843102225 + } + ], + "validators": { + "range": { + "min": -1.0, + "max": 1.0 + } + } + }, "imm_rates": { "title": "Immigration rates over the time path.", "description": "Immigration rates over the time path.", @@ -4598,6 +4615,45 @@ } } }, + "imm_rates_preTP": { + "title": "Immigration rates in period before model start year", + "description": "Immigration rates in period before model start year.", + "section_1": "Demographic Parameters", + "notes": "", + "type": "float", + "number_dims": 1, + "value": [ + { + "value": + [6.21124324e-03, 6.69539127e-03, 7.24160599e-03, 7.59057707e-03, + 7.64614775e-03, 7.51951110e-03, 7.06681489e-03, 6.56429771e-03, + 6.07319754e-03, 5.50700033e-03, 4.94732741e-03, 4.54958259e-03, + 4.23425001e-03, 3.97403995e-03, 3.67752183e-03, 3.38021287e-03, + 3.00514338e-03, 2.71916344e-03, 2.41686669e-03, 2.13548533e-03, + 1.86850952e-03, 1.63930855e-03, 1.45713458e-03, 1.34212247e-03, + 1.27164811e-03, 9.88138261e-04, 8.42267890e-04, 7.97801211e-04, + 6.63231805e-04, 4.87948389e-04, 4.32979019e-04, 5.60040782e-04, + 4.55369916e-04, 4.13717216e-04, 4.55029226e-04, 4.52104109e-04, + 5.58645176e-04, 6.14631276e-04, 5.66684121e-04, 5.11408015e-04, + 4.68968748e-04, 4.86541731e-04, 1.87084644e-04, -2.33525940e-05, + 1.21315737e-03, 1.05324177e-03, 8.29011129e-07, -1.80456466e-04, + 6.44696838e-04, 1.08838720e-03, 6.18137660e-04, 4.41473826e-04, + 6.10186803e-04, 6.58633423e-04, 5.63207404e-04, 5.34295717e-04, + 6.00114993e-04, 6.29050763e-04, 8.43098879e-04, 9.20493620e-04, + 9.68908696e-04, 1.26374022e-03, 1.60616956e-03, 1.72840304e-03, + 1.50429547e-03, 1.03459297e-03, 8.01449649e-04, 9.28625040e-04, + 1.60993336e-03, 2.23227102e-03, 3.53451125e-03, 5.32480275e-03, + 7.84622608e-03, 1.08220386e-02, 1.30810434e-02, 1.64519720e-02, + 1.88010017e-02, 2.12725905e-02, 2.17955114e-02, 1.92400972e-02] + } + ], + "validators": { + "range": { + "min": -1.0, + "max": 1.0 + } + } + }, "rho": { "title": "Age-specific mortality rates.", "description": "Age-specific mortality rates.", @@ -4617,6 +4673,44 @@ } } }, + "rho_preTP": { + "title": "Age-specific mortality rates.", + "description": "Age-specific mortality rates.", + "section_1": "Demographic Parameters", + "notes": "", + "type": "float", + "number_dims": 1, + "value": [ + { + "value": [7.32853396e-04, 8.17119742e-04, 8.78872967e-04, 9.11652521e-04, + 9.22980233e-04, 9.27340681e-04, 9.35710239e-04, 9.50095370e-04, + 9.75973557e-04, 1.01035546e-03, 1.04974597e-03, 1.08965056e-03, + 1.12905646e-03, 1.16544560e-03, 1.20432393e-03, 1.25117739e-03, + 1.31100453e-03, 1.38282557e-03, 1.46963538e-03, 1.57292023e-03, + 1.69015991e-03, 1.82531513e-03, 1.98785868e-03, 2.18224842e-03, + 2.40545835e-03, 2.64949570e-03, 2.91235990e-03, 3.19652189e-03, + 3.50144041e-03, 3.82658531e-03, 4.17597720e-03, 4.54600008e-03, + 4.92938669e-03, 5.32344068e-03, 5.73216261e-03, 6.18294920e-03, + 6.66506301e-03, 7.14289275e-03, 7.60533472e-03, 8.07929361e-03, + 8.60718380e-03, 9.21954043e-03, 9.91697732e-03, 1.07137947e-02, + 1.16186807e-02, 1.26511756e-02, 1.38062007e-02, 1.50692314e-02, + 1.64399087e-02, 1.79462854e-02, 1.96828400e-02, 2.16430913e-02, + 2.37466335e-02, 2.59837193e-02, 2.84282734e-02, 3.12585980e-02, + 3.45111290e-02, 3.80870625e-02, 4.19924732e-02, 4.63268386e-02, + 5.13310180e-02, 5.70418578e-02, 6.33201548e-02, 7.01769632e-02, + 7.77688358e-02, 8.62903314e-02, 9.59061899e-02, 1.06724255e-01, + 1.18798575e-01, 1.32125256e-01, 1.46682823e-01, 1.62432478e-01, + 1.79326966e-01, 1.97313414e-01, 2.16339448e-01, 2.35357633e-01, + 2.54059844e-01, 2.72119171e-01, 2.89195180e-01, 1.00000000e+00] + } + ], + "validators": { + "range": { + "min": 0.0, + "max": 1.0 + } + } + }, "etr_params": { "title": "Effective tax rate function parameters.", "description": "Effective tax rate function parameters.", diff --git a/ogcore/demographics.py b/ogcore/demographics.py index ce8a3e660..196f9bd93 100644 --- a/ogcore/demographics.py +++ b/ogcore/demographics.py @@ -339,7 +339,6 @@ def get_pop( infmort_rates=None, imm_rates=None, initial_pop=None, - pre_pop_dist=None, country_id=UN_COUNTRY_CODE, start_year=START_YEAR, end_year=END_YEAR, @@ -368,8 +367,6 @@ def get_pop( and model age initial_pop_data (Pandas DataFrame): initial population data for the first year of model calibration (start_year) - pre_pop_dist (Numpy array): population distribution for the year - before the initial year for calibration country_id (str): country id for UN data start_year (int): start year data end_year (int): end year for data @@ -377,33 +374,11 @@ def get_pop( Returns: pop_2D (Numpy array): population distribution over T0 periods - pre_pop (Numpy array): population distribution one year before - initial year for calibration of omega_S_preTP """ # Generate time path of the nonstationary population distribution # Get path up to end of data year pop_2D = np.zeros((end_year + 2 - start_year, E + S)) if infer_pop: - if pre_pop_dist is None: - pre_pop_data = get_un_data( - "47", - country_id=country_id, - start_year=start_year - 1, - end_year=start_year - 1, - ) - if download_path: - pre_pop_data.to_csv( - os.path.join(download_path, "raw_pre_pop_data_UN.csv"), - index=False, - ) - pre_pop_sample = pre_pop_data[ - (pre_pop_data["age"] >= min_age) - & (pre_pop_data["age"] <= max_age) - ] - pre_pop = pre_pop_sample.value.values - pre_pop_dist = pop_rebin(pre_pop, E + S) - else: - pre_pop = pre_pop_dist if initial_pop is None: initial_pop_data = get_un_data( "47", @@ -458,38 +433,17 @@ def get_pop( pop = pop_data_sample.value.values # Generate the current population distribution given that E+S might # be less than max_age-min_age+1 - # age_per_EpS = np.arange(1, E + S + 1) pop_EpS = pop_rebin(pop, E + S) pop_2D[y - start_year, :] = pop_EpS - # get population distribution one year before initial year for - # calibration of omega_S_preTP - pre_pop_data = get_un_data( - "47", - country_id=country_id, - start_year=start_year - 1, - end_year=start_year - 1, - ) - pre_pop_sample = pre_pop_data[ - (pre_pop_data["age"] >= min_age) & (pre_pop_data["age"] <= max_age) - ] - pre_pop = pre_pop_sample.value.values - if download_path: np.savetxt( os.path.join(download_path, "population_distribution.csv"), pop_2D, delimiter=",", ) - np.savetxt( - os.path.join( - download_path, "pre_period_population_distribution.csv" - ), - pre_pop, - delimiter=",", - ) - return pop_2D, pre_pop + return pop_2D def pop_rebin(curr_pop_dist, totpers_new): @@ -715,7 +669,6 @@ def get_pop_objs( imm_rates=None, infer_pop=False, pop_dist=None, - pre_pop_dist=None, country_id=UN_COUNTRY_CODE, initial_data_year=START_YEAR - 1, final_data_year=START_YEAR + 2, @@ -746,9 +699,6 @@ def get_pop_objs( infer_pop (bool): =True if want to infer the population pop_dist (array_like): user provided population distribution, dimensions are T0+1 x E+S - pre_pop_dist (array_like): user provided population distribution - for the year before the initial year for calibration, - length E+S country_id (str): country id for UN data initial_data_year (int): initial year of data to use (not relevant if have user provided data) @@ -775,9 +725,12 @@ def get_pop_objs( path, length T + S """ + start_data_year = initial_data_year - 1 # grab data from one year + T = T + 1 # add one period to T to account for period -1 pop + # before initial so have pre-start year population distribution # TODO: this function does not generalize with T. # It assumes one model period is equal to one calendar year in the - # time dimesion (it does adjust for S, however) + # time dimension (it does adjust for S, however) T0 = ( final_data_year - initial_data_year + 1 ) # number of periods until constant fertility and mortality rates @@ -788,8 +741,8 @@ def get_pop_objs( final_data_year, ) assert E + S <= max_age - min_age + 1 - assert initial_data_year >= 2011 and initial_data_year <= 2100 - 1 - assert final_data_year >= 2011 and final_data_year <= 2100 - 1 + assert initial_data_year >= 2012 and initial_data_year <= 2100 - 1 + assert final_data_year >= 2012 and final_data_year <= 2100 - 1 # Ensure that the last year of data used is before SS transition assumed # Really, it will need to be well before this assert final_data_year > initial_data_year @@ -810,7 +763,7 @@ def get_pop_objs( min_age, max_age, country_id, - initial_data_year, + start_data_year, final_data_year, download_path=download_path, ) @@ -839,7 +792,7 @@ def get_pop_objs( min_age, max_age, country_id, - initial_data_year, + start_data_year, final_data_year, download_path=download_path, ) @@ -877,7 +830,7 @@ def get_pop_objs( initial_pop = pop_dist[0, :].reshape(1, pop_dist.shape[-1]) else: initial_pop = None - pop_2D, pre_pop = get_pop( + pop_2D = get_pop( E, S, min_age, @@ -888,40 +841,32 @@ def get_pop_objs( infmort_rates, imm_rates, initial_pop, - pre_pop_dist, country_id, - initial_data_year, + start_data_year, final_data_year, download_path=download_path, ) else: - pop_2D, pre_pop = get_pop( + pop_2D = get_pop( E, S, min_age, max_age, country_id=country_id, - start_year=initial_data_year, + start_year=start_data_year, end_year=final_data_year, download_path=download_path, ) else: # Check first dims of pop_dist as input by user - print("T0 = ", T0) assert pop_dist.shape[0] == T0 + 1 # population needs to be # one year longer in order to find immigration rates assert pop_dist.shape[-1] == E + S - # Check that pre_pop specified - assert pre_pop_dist is not None - assert pre_pop_dist.shape[0] == pop_dist.shape[1] - pre_pop = pre_pop_dist # Create 2D array of population distribution pop_2D = np.zeros((T0 + 1, E + S)) for t in range(T0 + 1): pop_EpS = pop_rebin(pop_dist[t, :], E + S) pop_2D[t, :] = pop_EpS - # Get percentage distribution for S periods for pre-TP period - pre_pop_EpS = pop_rebin(pre_pop, E + S) # Get immigration rates if not provided if imm_rates is None: imm_rates_orig = get_imm_rates( @@ -933,7 +878,7 @@ def get_pop_objs( infmort_rates, pop_2D, country_id, - initial_data_year, + start_data_year, final_data_year, download_path=download_path, ) @@ -957,26 +902,13 @@ def get_pop_objs( ) # If the population distribution was given, check it for consistency # with the fertility, mortality, and immigration rates - # if pop_dist is not None and not infer_pop: - # len_pop_dist = pop_dist.shape[0] - # pop_counter_2D = np.zeros((len_pop_dist, E + S)) len_pop_dist = pop_2D.shape[0] pop_counter_2D = np.zeros((len_pop_dist, E + S)) # set initial population distribution in the counterfactual to # the first year of the user provided distribution - # pop_counter_2D[0, :] = pop_dist[0, :] pop_counter_2D[0, :] = pop_2D[0, :] for t in range(1, len_pop_dist): # find newborns next period - # newborns = np.dot(fert_rates[t - 1, :], pop_counter_2D[t - 1, :]) - - # pop_counter_2D[t, 0] = ( - # 1 - infmort_rates[t - 1] - # ) * newborns + imm_rates[t - 1, 0] * pop_counter_2D[t - 1, 0] - # pop_counter_2D[t, 1:] = ( - # pop_counter_2D[t - 1, :-1] * (1 - mort_rates[t - 1, :-1]) - # + pop_counter_2D[t - 1, 1:] * imm_rates_orig[t - 1, 1:] - # ) newborns = np.dot(fert_rates[t - 1, :], pop_counter_2D[t - 1, :]) pop_counter_2D[t, 0] = ( @@ -987,98 +919,8 @@ def get_pop_objs( + pop_counter_2D[t - 1, 1:] * imm_rates_orig[t - 1, 1:] ) # Check that counterfactual pop dist is close to pop dist given - # assert np.allclose(pop_counter_2D, pop_dist) assert np.allclose(pop_counter_2D, pop_2D) - # """" - # CHANGE - in OG-Core, we are implicitly assuming pre-TP rates of mortality, - # fertility, and immigration are the same as the period 0 rates. - - # So let's just infer the pre-pop_dist from those. - # """ - # pop1 = pop_2D[0, :] - # fert0 = fert_rates[0, :] - # mort0 = mort_rates[0, :] - # infmort0 = infmort_rates[0] - # imm0 = imm_rates_orig[0, :] - # pre_pop_guess = pop1.copy() - - # # I can't solve this analytically, so set up a system of equation - # # to solve - # def pre_pop_solve(pre_pop_guess, pop1, fert0, mort0, infmort0, imm0): - # pre_pop = pre_pop_guess - # errors = np.zeros(E + S) - # errors[0] = pop1[0] - ( - # (1 - infmort0) * (fert0 * pre_pop).sum() + imm0[0] * pre_pop[0] - # ) - # errors[1:] = pop1[1:] - ( - # pre_pop[:-1] * (1 - mort0[:-1]) + pre_pop[1:] * imm0[1:] - # ) - # # print("Max error = ", np.abs(errors).max()) - # return errors - - # opt_res = opt.root( - # pre_pop_solve, - # pre_pop_guess, - # args=(pop1, fert0, mort0, infmort0, imm0), - # method="lm", - # ) - # pre_pop = opt_res.x - # print( - # "Success? ", - # opt_res.success, - # ", Max diff = ", - # np.abs(opt_res.fun).max(), - # ) - # pre_pop_EpS = pop_rebin(pre_pop, E + S) - - # # Check result - # initial_pop_counter = np.zeros(E + S) - # newborns = (fert_rates[0, :] * pre_pop[:]).sum() - # initial_pop_counter[0] = ( - # 1 - infmort_rates[0] - # ) * newborns + imm_rates_orig[0, 0] * pre_pop[0] - # initial_pop_counter[1:] = ( - # pre_pop[:-1] * (1 - mort_rates[0, :-1]) - # + pre_pop[1:] * imm_rates_orig[0, 1:] - # ) - # # Test that using pre pop get to pop in period 1 - # print("Max diff = ", np.abs(pop_2D[0, :] - initial_pop_counter).max()) - # # assert np.allclose(initial_pop_counter, pop_2D[0, :]) - - """ - NEW CODE - use the actual UN historical data instead of solving backwards - """ - # Get percentage distribution for S periods for pre-TP period - # Use the actual UN historical data instead of solving backwards - # pre_pop_EpS = pop_rebin(pre_pop, E + S) -- this assignment is - # above on line 924, but keep for clarity - - # Check result - # Verify that the UN pre-period data is reasonably consistent - # with the period 0 population using the demographic transition equations - initial_pop_counter = np.zeros(E + S) - newborns = (fert_rates[0, :] * pre_pop_EpS[:]).sum() - initial_pop_counter[0] = ( - 1 - infmort_rates[0] - ) * newborns + imm_rates_orig[0, 0] * pre_pop_EpS[0] - initial_pop_counter[1:] = ( - pre_pop_EpS[:-1] * (1 - mort_rates[0, :-1]) - + pre_pop_EpS[1:] * imm_rates_orig[0, 1:] - ) - - max_diff = np.abs(pop_2D[0, :] - initial_pop_counter).max() - print("Pre-period population verification: Max diff = ", max_diff) - - if max_diff > 100_000: - print( - "WARNING: Large difference between UN pre-period population " - + "and period 0 population ({:.2f}). ".format(max_diff) - + "This may indicate inconsistencies in the data or " - + "immigration rate calculations, but using UN historical " - + "data as it is more reliable than backward-solved estimates." - ) - # Create the transition matrix for the population distribution # from T0 going forward (i.e., past when we have data on forecasts) OMEGA_orig = np.zeros((E + S, E + S)) @@ -1112,7 +954,8 @@ def get_pop_objs( # steady-state distribution by adjusting immigration rates, holding # constant mortality, fertility, and SS growth rates imm_tol = 1e-14 - fixper = int(1.5 * S + T0) + fixper = int(1.5 * S) + assert fixper > T0 # ensure that we are fixing period after data omega_SSfx = omega_path_lev[fixper, :] / omega_path_lev[fixper, :].sum() imm_objs = ( fert_rates[fixper, :], @@ -1137,15 +980,11 @@ def get_pop_objs( omega_path_S[fixper, :].reshape((1, S)), (T + S - fixper, 1) ) g_n_path = np.zeros(T + S) - g_n_path[1:] = ( + g_n_path[:-1] = ( omega_path_lev[1:, -S:].sum(axis=1) - omega_path_lev[:-1, -S:].sum(axis=1) ) / omega_path_lev[:-1, -S:].sum(axis=1) - g_n_path[0] = ( - omega_path_lev[0, -S:].sum() - pre_pop_EpS[-S:].sum() - ) / pre_pop_EpS[-S:].sum() g_n_path[fixper + 1 :] = g_n_SS - omega_S_preTP = pre_pop_EpS[-S:] / pre_pop_EpS[-S:].sum() imm_rates_mat = np.concatenate( ( imm_rates_orig[:fixper, E:], @@ -1304,13 +1143,16 @@ def get_pop_objs( # Return objects in a dictionary pop_dict = { - "omega": omega_path_S, + "omega": omega_path_S[1:, :], "g_n_ss": g_n_SS, "omega_SS": omega_SSfx[-S:] / omega_SSfx[-S:].sum(), - "rho": mort_rates_S, - "g_n": g_n_path, - "imm_rates": imm_rates_mat, - "omega_S_preTP": omega_S_preTP, + "rho": mort_rates_S[1:, :], + "g_n": g_n_path[1:], + "imm_rates": imm_rates_mat[1:, :], + "omega_S_preTP": omega_path_S[0, :], + "imm_rates_preTP": imm_rates_mat[0, :], + "rho_preTP": mort_rates_S[0, :], + "g_n_preTP": g_n_path[0], } return pop_dict diff --git a/ogcore/parameters.py b/ogcore/parameters.py index 0134d21fc..a3d39a76c 100644 --- a/ogcore/parameters.py +++ b/ogcore/parameters.py @@ -317,6 +317,7 @@ def compute_default_params(self): # for constant demographics if self.constant_demographics: self.g_n_ss = 0.0 + self.g_n_preTP = 0.0 self.g_n = np.zeros(self.T + self.S) surv_rate = np.ones_like(self.rho) - self.rho surv_rate1 = np.ones((self.S,)) # prob start at age S @@ -329,6 +330,8 @@ def compute_default_params(self): np.reshape(self.omega_SS, (1, self.S)), (self.T + self.S, 1) ) self.omega_S_preTP = self.omega_SS + self.imm_rates_preTP = np.zeros(self.S) + self.rho_preTP = self.rho[0, :] # Create time series of stationarized UBI transfers self.ubi_nom_array = self.get_ubi_nom_objs() diff --git a/tests/test_aggregates.py b/tests/test_aggregates.py index 9cfdf8cb8..2cae0ce23 100644 --- a/tests/test_aggregates.py +++ b/tests/test_aggregates.py @@ -151,6 +151,10 @@ def test_get_I(b_splus1, K_p1, K, p, method, expected): "omega": np.ones((160, 40)) / 40, "omega_SS": np.ones(40) / 40, "imm_rates": np.zeros((160, 40)), + "omega_S_preTP": np.ones(40) / 40, + "imm_rates_preTP": np.zeros(40), + "rho_preTP": rho_vec[0, :], + "g_n_preTP": 0.01, } # update parameters instance with new values for test p.update_specifications(new_param_values) @@ -173,7 +177,7 @@ def test_get_I(b_splus1, K_p1, K, p, method, expected): expected2 = B_test.sum(1).sum(1) / ( 1.0 + np.hstack((p.g_n[1 : p.T], p.g_n_ss)) ) -expected3 = B_test[0, :, :].sum() / (1.0 + p.g_n[0]) +expected3 = B_test[0, :, :].sum() / (1.0 + p.g_n_preTP) test_data = [ (b[-1, :, :], p, "SS", False, expected1), (b, p, "TPI", False, expected2), @@ -210,7 +214,10 @@ def test_get_B(b, p, method, PreTP, expected): "omega": np.ones((160, 40)) / 40, "omega_SS": np.ones(40) / 40, "imm_rates": np.zeros((160, 40)), - "rho": rho_vec.tolist(), + "omega_S_preTP": np.ones(40) / 40, + "imm_rates_preTP": np.zeros(40), + "rho_preTP": rho_vec[0, :], + "g_n_preTP": 0.01, } # update parameters instance with new values for test p.update_specifications(new_param_values) @@ -223,6 +230,7 @@ def test_get_B(b, p, method, PreTP, expected): np.reshape(p.rho[0, :] * pop, (p.T, p.S, 1)), (1, 1, p.J) ) growth_adj = (1.0 + r) / (1.0 + p.g_n[: p.T]) +growth_adj_preTP = (1.0 + r[0]) / (1.0 + p.g_n_preTP) expected1 = BQ_presum[-1, :, :].sum(0) * growth_adj[-1] expected2 = BQ_presum[-1, :, 1].sum(0) * growth_adj[-1] @@ -230,8 +238,8 @@ def test_get_B(b, p, method, PreTP, expected): np.reshape(growth_adj, (p.T, 1)), (1, p.J) ) expected4 = BQ_presum[:, :, 1].sum(1) * growth_adj -expected5 = BQ_presum[0, :, :].sum(0) * growth_adj[0] -expected6 = BQ_presum[0, :, 1].sum(0) * growth_adj[0] +expected5 = BQ_presum[0, :, :].sum(0) * growth_adj_preTP +expected6 = BQ_presum[0, :, 1].sum(0) * growth_adj_preTP p2 = copy.deepcopy(p) p2.use_zeta = True diff --git a/tests/test_demographics.py b/tests/test_demographics.py index 7ea0e2167..48f6d1fde 100644 --- a/tests/test_demographics.py +++ b/tests/test_demographics.py @@ -24,10 +24,6 @@ pop_dist = np.loadtxt( os.path.join(data_dir, "population_distribution.csv"), delimiter="," ) -pre_pop_dist = np.loadtxt( - os.path.join(data_dir, "pre_period_population_distribution.csv"), - delimiter=",", -) @pytest.mark.local @@ -77,7 +73,6 @@ def test_get_pop_objs(): imm_rates=imm_rates, infer_pop=True, pop_dist=pop_dist[0, :].reshape(1, E + S), - pre_pop_dist=pre_pop_dist, initial_data_year=start_year - 1, final_data_year=start_year, GraphDiag=False, @@ -107,7 +102,6 @@ def test_pop_smooth(): imm_rates=imm_rates, infer_pop=True, pop_dist=pop_dist[0, :].reshape(1, E + S), - pre_pop_dist=pre_pop_dist, initial_data_year=start_year - 1, final_data_year=start_year, country_id="840", @@ -155,7 +149,6 @@ def test_pop_growth_smooth(): imm_rates=imm_rates, infer_pop=True, pop_dist=pop_dist[0, :].reshape(1, E + S), - pre_pop_dist=pre_pop_dist, initial_data_year=start_year - 1, final_data_year=start_year, country_id="840", @@ -363,7 +356,6 @@ def test_custom_series_fail(): "47", start_year=start_year - 1, end_year=start_year - 1 ) pop = df[(df.age < 100) & (df.age >= 0)].value.values - pre_pop_dist = demographics.pop_rebin(pop, E + S) pop_dict = demographics.get_pop_objs( E, S, @@ -375,7 +367,6 @@ def test_custom_series_fail(): infmort_rates=infmort_rates, imm_rates=imm_rates, pop_dist=pop_dist, - pre_pop_dist=pre_pop_dist, initial_data_year=start_year, final_data_year=start_year + 1, GraphDiag=False, @@ -405,7 +396,6 @@ def test_SS_dist(): imm_rates=imm_rates, infer_pop=True, pop_dist=pop_dist[0, :].reshape(1, E + S), - pre_pop_dist=pre_pop_dist, initial_data_year=start_year - 1, final_data_year=start_year, GraphDiag=False, @@ -438,7 +428,6 @@ def test_time_path_length(): imm_rates=imm_rates, infer_pop=True, pop_dist=pop_dist[0, :].reshape(1, E + S), - pre_pop_dist=pre_pop_dist, initial_data_year=start_year - 1, final_data_year=start_year, GraphDiag=False, @@ -450,7 +439,7 @@ def test_time_path_length(): assert pop_dict["rho"].shape[0] == T + S -# test of get pop when infer population, but don't pass initial pop or pre_pop +# test of get pop when infer population, but don't pass initial pop @pytest.mark.local def test_infer_pop_nones(): """ @@ -501,7 +490,6 @@ def test_infer_pop_nones(): imm_rates=imm_rates, infer_pop=True, pop_dist=None, - pre_pop_dist=None, initial_data_year=start_year, final_data_year=start_year + 1, GraphDiag=False, @@ -548,29 +536,36 @@ def test_data_download(tmpdir): pop_dist = np.loadtxt( os.path.join(tmpdir, "population_distribution.csv"), delimiter="," ) - pre_pop_dist = np.loadtxt( - os.path.join(tmpdir, "pre_period_population_distribution.csv"), - delimiter=",", - ) pop_dict2 = demographics.get_pop_objs( E, S, T, 0, 99, - fert_rates=fert_rates, - mort_rates=mort_rates, - infmort_rates=infmort_rates, - imm_rates=imm_rates, + fert_rates=fert_rates[:, :], + mort_rates=mort_rates[:, :], + infmort_rates=infmort_rates[:], + imm_rates=imm_rates[:, :], infer_pop=True, pop_dist=pop_dist[0, :].reshape(1, E + S), - pre_pop_dist=pre_pop_dist, initial_data_year=start_year, - final_data_year=start_year + 1, + final_data_year=start_year + 2, GraphDiag=False, ) # Assert that the two pop_dicts are the same for key in pop_dict: print(key) - assert np.allclose(pop_dict[key], pop_dict2[key]) + print("Diff =", np.abs(pop_dict[key] - pop_dict2[key]).max()) + if key == "imm_rates": + # print the max diff for each T + for t in range(pop_dict[key].shape[0]): + print( + "Max diff for imm_rates at T = ", + t, + " is ", + np.abs(pop_dict[key][t, :] - pop_dict2[key][t, :]).max(), + ) + for key in pop_dict: + print(key) + assert np.allclose(pop_dict[key], pop_dict2[key], atol=7e-5) diff --git a/tests/testing_params.json b/tests/testing_params.json index d38869512..2f1206802 100644 --- a/tests/testing_params.json +++ b/tests/testing_params.json @@ -2334,5 +2334,22 @@ 0.03104479, 0.02942136, 0.02786169, 0.02681268, 0.02060319, 0.0197469 , 0.01730843, 0.01474155, 0.01285935, 0.01119112, 0.00984811, 0.00801472, 0.00696934, 0.00592803, 0.00454995, - 0.00352386, 0.00256576, 0.00167457, 0.00106337, 0.00061474] + 0.00352386, 0.00256576, 0.00167457, 0.00106337, 0.00061474], +"rho_preTP": [0.0025781 , 0.00268902, 0.00282206, 0.00303204, 0.00333702, + 0.00373181, 0.00435139, 0.00521627, 0.00608546, 0.00715621, + 0.0084698 , 0.00975143, 0.01134659, 0.01332463, 0.01493943, + 0.01647184, 0.0181378 , 0.01991245, 0.02238255, 0.02570651, + 0.02896702, 0.03292571, 0.0381694 , 0.04358988, 0.05031253, + 0.05940869, 0.06896693, 0.08087685, 0.09692589, 0.11366492, + 0.13429498, 0.16117335, 0.1886646 , 0.22230383, 0.26433221, + 0.30426546, 0.35003558, 0.39994367, 0.43838896, 1.0], +"imm_rates_preTP": [0.00172906, 0.00116727, -0.00016137, -0.01018377, -0.01893683, -0.00107257 +, 0.00197262, 0.01956533, 0.01741252, -0.00494443, -0.00777153, 0.01939821 +, 0.01978244, 0.00747581, 0.00741143, 0.00150033, -0.00766893, -0.0058054 +, -0.00994736, -0.02296918, -0.00731107, -0.00217819, -0.02644512, -0.04922738 +, -0.01088844, -0.03118354, -0.03035488, -0.00925152, -0.02137357, 0.00121973 +, 0.00414307, -0.00055469, 0.02895412, 0.03014472, 0.00587107, 0.05014318 +, 0.05761331, 0.01828307, 0.08218118, 0.06926356], +"g_n_preTP": 7.96404268e-03 } +