Problems with the replication of Castaneda et al. (2003) on CPU

1

Good afternoon

First of all, thanks @robertdkirkby for having developed the toolkit and having made it publicly available!

I am writing here because I’ve had some problems working with the toolkit. For an undergraduate project, I’m trying to replicate the model used in Castaneda et al. (2003) and calibrate it on the Italian economy. I, therefore, found very useful the replication of the paper you have made available. However, I don’t have access to a machine that can run MATLAB with a GPU. I changed the options in order to make the code compatible with a CPU and eliminated any instance in which an array would have been stored on GPU.

When I try to run the code, however, I get an error. This happens both for my modified version for Italy and the original replication code. Below I have copied the error message I receive. Moreover, in this GoogleDrive folder, there are the codes I’m using for both the Italian model and for running Castaneda on CPU. Do you know what could be wrong, and how I could get the model to work?

Thank you in advance!

Andrea

https://drive.google.com/drive/folders/10_NXMBkSV3nnozl9gOFt8xNN7A1Yl0Zz?usp=sharing

Not enough input arguments.
Error in
CPU_Italy_Model>@(d1_val,d2_val,a_val,z_val,r,sigma1,sigma2,chi,elle,theta,delta,e1,e2,e3,e4,e5,omega,a0,a1,a2,a3)CPU_Italy_Model_ReturnFn(d1_val,d2_val,a_val,z_val,r,sigma1,sigma2,chi,elle,theta,delta,e1,e2,e3,e4,e5,omega,a0,a1,a2,a3)
(line 194)
ReturnFn=@(d1_val, d2_val, a_val, z_val,r,sigma1,sigma2,chi,elle,theta,delta,e1,e2,e3,e4,e5,omega,a0,a1,a2,a3)
CPU_Italy_Model_ReturnFn(d1_val, d2_val, a_val, z_val,r,sigma1,sigma2,chi,elle,theta,delta,e1,e2,e3,e4,e5,omega,a0,a1,a2,a3);
Error in CreateReturnFnMatrix_Case2_Disc (line 58)

  •            Fmatrix(i1,i2,i3)=ReturnFn(d_gridvals(i1,:),a_gridvals(i2,:),z_gridvals(i3,:));*

Error in ValueFnIter_Case2 (line 126)

  •    ReturnMatrix=CreateReturnFnMatrix_Case2_Disc(ReturnFn, n_d, n_a, n_z, d_grid, a_grid, z_grid, vfoptions.parallel);*

Error in HeteroAgentStationaryEqm_Case2_subfn (line 12)
[~, Policy]=ValueFnIter_Case2(V0Kron, n_d, n_a, n_s, d_grid, a_grid, s_grid, pi_s, Phi_aprimeKron, Case2_Type, ReturnFn, Parameters,
DiscountFactorParamNames, ReturnFnParamNames, PhiaprimeParamNames, vfoptions);
Error in
HeteroAgentStationaryEqm_Case2>@§HeteroAgentStationaryEqm_Case2_subfn(p,V0Kron,n_d,n_a,n_s,pi_s,d_grid,a_grid,s_grid,Phi_aprimeKron,Case2_Type,ReturnFn,FnsToEvaluate,GeneralEqmEqns,Parameters,DiscountFactorParamNames,ReturnFnParamNames,PhiaprimeParamNames,FnsToEvaluateParamNames,GeneralEqmEqnParamNames,GEPriceParamNames,heteroagentoptions,simoptions,vfoptions)
(line 69)
GeneralEqmConditionsFn=@§ HeteroAgentStationaryEqm_Case2_subfn(p, V0Kron, n_d, n_a, n_s, pi_s, d_grid, a_grid, s_grid, Phi_aprimeKron,
Case2_Type, ReturnFn, FnsToEvaluate, GeneralEqmEqns, Parameters, DiscountFactorParamNames, ReturnFnParamNames, PhiaprimeParamNames,
FnsToEvaluateParamNames, GeneralEqmEqnParamNames, GEPriceParamNames,heteroagentoptions, simoptions, vfoptions)
Error in fminsearch (line 200)
fv(:,1) = funfcn(x,varargin{:});
Error in HeteroAgentStationaryEqm_Case2 (line 80)

  • [p_eqm,GeneralEqmConditions]=fminsearch(GeneralEqmConditionsFn,p0);*
    Error in CPU_Italy_Model (line 250)
  • [p_eqm,p_eqm_index,MarketClearance]=HeteroAgentStationaryEqm_Case2(V0, n_d, n_a, n_z, 0, pi_z, d_grid, a_grid, z_grid,Phi_aprimeMatrix,*
  • Case2_Type, ReturnFn, FnsToEvaluate, GeneralEqmEqns, Params, DiscountFactorParamNames, ReturnFnParamNames, PhiaprimeParamNames,*
  • FnsToEvaluateParamNames, GeneralEqmEqnParamNames, GEPriceParamNames,heteroagentoptions, simoptions, vfoptions);*
    Error in run (line 91)
    evalin(‘caller’, strcat(script, ‘;’));
2 Likes
2

Dear Andrea,

was my fault. I rarely use cpu so had let the Case2 ValueFnIter command get rusty with cpu. Have pushed some changes to github that fix this. At least when I ran your Castaneda cpu codes (did not try the Italy ones).

I have done a lot of double and triple checking the outputs when using gpu, but not with cpu. So I suggest that you first run using cpu and check the outputs against those from gpu (you can find them in the relevant appendix of an (in progress) paper I have on replication) to double check that cpu is doing what I should be (I am confident that it is, just that I have not done any serious testing of it). Then do your version using Italy.
[“Quantitative Macro: Lessons Learnt from Eleven Replications”, link from http://www.robertdkirkby.com/research/ ]

Just a warning it will be a lot slower on cpu; hence why I concentrate on developing on the gpu. Also there is no reason not to use parallel=1 (parallel cpus, essentially all computers have them nowadays) rather than 0 (single cpu). That will be faster.

Let me know if anything else doesn’t work :slight_smile:

Yours Sincerely
Robert Kirkby

3 Likes
3

Also, when switching to Italy, there are a few lines around line 200 that are currently commented out that do the value function and agent distribution without calculating the general equilibrium. You might find it useful to look at how your recalibration is changing those before you calculate the whole general eqm.

1 Like
4

Thank you very much for the quick answer and for having dedicated time to this! I will check tomorrow morning (UK time) if it works and will let you know.

Yours,

Andrea

1 Like
5

Good afternoon again!

I tried to run the code, but unfortunately I was unsuccessful. When running the general equilibrium function for heterogenous agents, the code brought me to ind2sub_homemade_gpu.m. This function failed because it requires a GPU to store an array.

I tried to override the problem by preventing the array to be stored on GPU. When I run again the main code, it now gives the error code attached below. The same error takes place both when running the CPU version of Castaneda (2003) and when running the Italian model. Do you know what could be the source of the problem?

On another note, the value function iteration check around line 200 now works fine for both models, while it didn’t before!

Yours,

Andrea

Error code:

Index in position 1 exceeds array bounds.

Error in CreateGridvals_Policy (line 176)

  •        d_ind=PolicyIndexes(j1,j2);*

Error in EvalFnOnAgentDist_AggVars_Case2 (line 41)

  • [d_gridvals, ~]=CreateGridvals_Policy(PolicyIndexes,n_d,n_a,n_a,n_z,d_grid,a_grid,2, 2);*

Error in HeteroAgentStationaryEqm_Case2_subfn (line 17)
AggVars=EvalFnOnAgentDist_AggVars_Case2(StationaryDistKron, Policy, FnsToEvaluate, Parameters, FnsToEvaluateParamNames, n_d, n_a,
n_s, d_grid, a_grid, s_grid, simoptions.parallel);

Error in
HeteroAgentStationaryEqm_Case2>@§HeteroAgentStationaryEqm_Case2_subfn(p,V0Kron,n_d,n_a,n_s,pi_s,d_grid,a_grid,s_grid,Phi_aprimeKron,Case2_Type,ReturnFn,FnsToEvaluate,GeneralEqmEqns,Parameters,DiscountFactorParamNames,ReturnFnParamNames,PhiaprimeParamNames,FnsToEvaluateParamNames,GeneralEqmEqnParamNames,GEPriceParamNames,heteroagentoptions,simoptions,vfoptions)
(line 69)
GeneralEqmConditionsFn=@§ HeteroAgentStationaryEqm_Case2_subfn(p, V0Kron, n_d, n_a, n_s, pi_s, d_grid, a_grid, s_grid,
Phi_aprimeKron, Case2_Type, ReturnFn, FnsToEvaluate, GeneralEqmEqns, Parameters, DiscountFactorParamNames, ReturnFnParamNames,
PhiaprimeParamNames, FnsToEvaluateParamNames, GeneralEqmEqnParamNames, GEPriceParamNames,heteroagentoptions, simoptions, vfoptions)

Error in fminsearch (line 200)
fv(:,1) = funfcn(x,varargin{:});

Error in HeteroAgentStationaryEqm_Case2 (line 80)

  • [p_eqm,GeneralEqmConditions]=fminsearch(GeneralEqmConditionsFn,p0);*

Error in CastanedaDiazGimenezRiosRull2003 (line 217)

  • [p_eqm,p_eqm_index,MarketClearance]=HeteroAgentStationaryEqm_Case2(V0, n_d, n_a, n_z, 0, pi_z, d_grid, a_grid,*
  • z_grid,Phi_aprimeMatrix, Case2_Type, ReturnFn, FnsToEvaluate, GeneralEqmEqns, Params, DiscountFactorParamNames,*
  • ReturnFnParamNames, PhiaprimeParamNames, FnsToEvaluateParamNames, GeneralEqmEqnParamNames,*
  • GEPriceParamNames,heteroagentoptions, simoptions, vfoptions);*
1 Like
6

Just pushed some updates to github that solve this (correcting CreateGridvals_Policy() to allow for just using cpu). That part should now work. Let me know if you come to a further error.

1 Like
7

Thank you!

It does work neatly now.

1 Like
9

Hello.

I tried to replicate your code for CDR 2003. What I first did was change the parallel=1 in all files.

One of the errors I had was the following:

Unrecognized function or variable ‘aprimePolicy_gridvals’.

Error in CreateGridvals_Policy (line 103)
aprimePolicy_gridvals=num2cell(aprimePolicy_gridvals);

Error in EvalFnOnAgentDist_AggVars_Case2 (line 72)
[d_gridvals, ~]=CreateGridvals_Policy(PolicyIndexes,n_d,n_a,n_a,n_z,d_grid,a_grid,2, 2);

Error in CastanedaDiazGimenezRiosRull2003 (line 136)
AggVars=EvalFnOnAgentDist_AggVars_Case2(StationaryDist, Policy, FnsToEvaluate, Params, [], n_d, n_a, n_z, d_grid, a_grid, z_grid, 2);

Could you please help me with this? Thank you very much!!

1 Like
10

I believe this should now be fixed. Please update VFI Toolkit and try running again. Apologies for that.

(I updated CreateGridvals_Policy(), an internal function you will never directly use, about 2 months ago, reworking to make it way faster. I forgot to test it on a Case 2 problem and so there was a ‘nothing’ error. Simple one line fix.)

1 Like
11

Hello, I uninstalled VFI Toolkit and installed it again from Matlab add-ons. I reran the code with the parallel modifications, and unfortunately, I got the same error after 1 hour and something of the process.

Thank you for any help or comment!

Pablo

1 Like
12

Hmm, Matlab (mathworks file exchange) only syncs to the github periodically. It says on the mathworks website that last sync was “Updated Thu, 01 Dec 2022”. I guess it will resync soon (it periodically automatically resyncs, I’m not sure how often)? I cannot find a way to force sync (did a quick search around but found nothing).

You can either wait for it to sync. Or directly download from github. Or you can just edit ‘CreateGrid_Policy’ to add what is now line 80 as:
aprimePolicy_gridvals=nan; % there is no aprime when using Case2
(as seen at:
https://github.com/vfitoolkit/VFIToolkit-matlab/blob/master/SubCodes/CreateGridvals_Policy.m
which is the fix that I made)

2 Likes
13

Thank you for the reply, Sir.
I did the change as indicated, the previous error did not appear and it ran for a couple of hours, but then I got this:

Starting parallel pool (parpool) using the ‘Processes’ profile …
Connected to the parallel pool (number of workers: 4).
If this number is not zero there is an problem with Phi_aprimeKron: 0
GeneralEqmConditionsFnOpt =

function_handle with value:

@(p)HeteroAgentStationaryEqm_Case2_subfn(p,n_d,n_a,n_s,l_p,pi_s,d_grid,a_grid,s_grid,Phi_aprimeKron,Case2_Type,ReturnFn,FnsToEvaluate,GeneralEqmEqns,Parameters,DiscountFactorParamNames,ReturnFnParamNames,PhiaprimeParamNames,FnsToEvaluateParamNames,GeneralEqmEqnParamNames,GEPriceParamNames,heteroagentoptions,simoptions,vfoptions)

Unrecognized field name “outputGEform”.

Error in HeteroAgentStationaryEqm_Case2_subfn (line 29)
if heteroagentoptions.outputGEform==0 % scalar

Error in HeteroAgentStationaryEqm_Case2>@(p)HeteroAgentStationaryEqm_Case2_subfn(p,n_d,n_a,n_s,l_p,pi_s,d_grid,a_grid,s_grid,Phi_aprimeKron,Case2_Type,ReturnFn,FnsToEvaluate,GeneralEqmEqns,Parameters,DiscountFactorParamNames,ReturnFnParamNames,PhiaprimeParamNames,FnsToEvaluateParamNames,GeneralEqmEqnParamNames,GEPriceParamNames,heteroagentoptions,simoptions,vfoptions) (line 103)
GeneralEqmConditionsFnOpt=@(p) HeteroAgentStationaryEqm_Case2_subfn(p, n_d, n_a, n_s,l_p, pi_s, d_grid, a_grid, s_grid, Phi_aprimeKron, Case2_Type, ReturnFn, FnsToEvaluate, GeneralEqmEqns, Parameters, DiscountFactorParamNames, ReturnFnParamNames, PhiaprimeParamNames, FnsToEvaluateParamNames, GeneralEqmEqnParamNames, GEPriceParamNames,heteroagentoptions, simoptions, vfoptions)

Error in fminsearch (line 201)
fv(:,1) = funfcn(x,varargin{:});

Error in HeteroAgentStationaryEqm_Case2 (line 116)
[p_eqm_vec,GeneralEqmConditions]=fminsearch(GeneralEqmConditionsFnOpt,p0);

Error in CastanedaDiazGimenezRiosRull2003 (line 145)
[p_eqm,p_eqm_index,MarketClearance]=HeteroAgentStationaryEqm_Case2(n_d, n_a, n_z, 0, pi_z, d_grid, a_grid, z_grid,Phi_aprimeMatrix, Case2_Type, ReturnFn, FnsToEvaluate, GeneralEqmEqns, Params, DiscountFactorParamNames, [], PhiaprimeParamNames, [], [], GEPriceParamNames,heteroagentoptions, simoptions, vfoptions);

I can wait for the update, but I’d like to replicate all the models from your paper.

Thank you very much.

Pablo

14

Have corrected this (is the file VFIToolkit-matlab/HeteroAgentStationaryEqm_Case2.m at master · vfitoolkit/VFIToolkit-matlab · GitHub )

If you just want to make sure things work I strongly suggest reducing the grid sizes (e.g. n_d=21 and n_a=301) and running that. On the cpu it is going to have a difficult time solving things with grids as large as those in the replication.

[the replications are deliberately using ‘overly’ large grids to ensure accuracy]

2 Likes
15

Update available; let’s see if it works now! Fingers crossed!

16

Hello. I am almost done replicating this CDR(2003) in a CPU. I did as you suggested. It takes some hours.
I got an error just after table 8:

Error using arrayfun
All of the input arguments must be of the same size and shape.
Previous inputs had size 8 in dimension 2. Input #4 has size 1

Error in EvalFnOnAgentDist_Grid_Case2 (line 208)
Values=arrayfun(FnToEvaluate, d1vals, d2vals, a1vals, z1vals, ParamCell{:});

Error in EvalFnOnAgentDist_RankTransitionProbabilities_Case2 (line 134)
Values=EvalFnOnAgentDist_Grid_Case2(FnsToEvaluate{jj}, FnToEvaluateParamsVec,PolicyValuesPermute,n_d,n_a,n_z,a_grid,z_grid,Parallel);

Any help on this? Please let me know if I arrived at the point I am spamming.

Thank you!

Pablo

17

Update: I also ran the original code in a GPU and got the same error.

18

I just ran this on gpu without error. [Not sure why you would have got an error. I double-checked that the versions on my harddrive are fully synced to the github repo.]

I am going to try rerunning on cpu, see if I get the same error as you there?

Apologies for taking so long to respond to this, both I and the computer were busy with other things :slight_smile:

I should hopefully get back to you about cpu next week.

[I understand what the error is saying is the issue, but I have all the inputs shaped in such a way that this error should not be thrown. arrayfun() has different functionality on the cpu and gpu, but in both cases the idea is that we evaluate the FnToEvaluate on the cross-product grid from the two decision variables (d1vals and d2vals) the endogenous state (a1vals) and the exogenous state (z1vals).]

1 Like
19

Thank you very much!
Just to be clear, when I ran it on a GPU provided by my school, the only change I made was the following:

n_l=21;
n_k=301;

I did this because when I tried to put the original values

n_l=51
n_k=1501

It ran out of memory.

Best!

Pablo

20

Changing the grid sizes should not cause an error (except of course out of memory ones when the grid is too big :wink: )

These replications all use grids that are probably ‘too big’, in the sense that I used deliberately excessive grids to make sure the replications were very accurate (I want the replications to be as accurate as I can get them).

1 Like
21

Probably, the issue on the cpu was that it is ‘forcing’ the gpu (which is obviously not possible).

On line 318 of CastanedaDiazGimenezRiosRull2003:

TransitionProbabilities=EvalFnOnAgentDist_RankTransitionProbabilities_Case2(t,NSims,StationaryDist, Policy,Phi_aprimeMatrix, Case2_Type, FnsToEvaluate3, Params,, n_d, n_a, n_z, d_grid, a_grid, z_grid, pi_z, 2, npoints);

The second last input (the 2) is ‘Parallel’ which when given the value of 2 is gpu parallelization. If you change it to 1 (parallel cpus) then that should mean it now works.

In the current VFI Toolkit the parallelization is typically automatically detected and set (and can be forced by the user using options). This replication was one of the early ones and so used a now outdated version of VFI Toolkit in which these things had to be declared manually.