Saturday, March 5, 2016

SIM2: Interactive spreadsheet implementing Godley & Lavoie's SIM model (now w/ sample period invariance! [Hint: the secret is expm(A*Ts) & its integral have Ts parameters!])

tl;dr version: Everything you need to know is right here (except when A is singular, in which case you'll need this too to find discrete time B) which tells you about converting back and forth from continuous to discrete time in a problem like this, namely, given the system:

dx/dt = A∙x + B∙u       

for vectors x(t) and u(t) and constant matrices A and B and u(t) held constant over each sample period T, we have the discrete time equivalent:
Notice that the discrete time constant matrix coefficients Ad and Bd are dependent on the sample period T.
end tl;dr

Success!! (I think). I changed my previous spreadsheet to make it truly invariant to the sample period (Ts in the green box below). Please refer to the previous version and the description underneath the spreadsheet for what the formulas are doing and the discrete time state space representation of the system, and to links to Jason and Ramanan's blogs (from whence this whole exercise originated). Below the spreadsheet on this page I describe what I'm doing differently. Feel free to play with the values: I only planned on people editing the green cells (in this embedded version), and if you get it irrevocably screwed up, just refresh the page (you can't hurt my original copy). Also you can download your own copy using the tool on the right side of the black toolbar running along the bottom of the sheet. Note that setting Ts < 0.1 will limit the data to less than 100 years (so the plot won't be filled out) because the number of time samples in my table is fixed at 1000.

Note the names and values with the modifier "alt" in them below don't affect anything: they're part of an experiment to get SIM2 to always satisfy all of G&L's accounting expressions. I've had success with the ΔH = G - T expression by using T-alt (in columns to the right, you have to scroll to see). I did add T-alt to the plot so you can see how it differs slightly from T. More in Update 5 below and here.

NOTE (update of 2016.04.22): in the above sheet to calculate the so-called experimental "alt" version of values use CTalt = (1-A)/Ts and DTalt = (Ts-B)/Ts and alpha2alt = CTalt/DTalt, where A and B are a function of the new sample period Ts (assuming the old one Ts0 = 1). These expressions come from T as a function of Ts := T(Ts) = T*Ts = G(Ts) - ΔH(Ts) = G[n]*Ts - H[n](Ts) + H[n-1] = G[n]*Ts - A*H[n-1] - B*G[n] + H[n-1] = (1-A)*H[n-1] + (Ts-B)*G[n] = CTalt*H[n-1] + DTalt*G[n], assuming H[n-1] didn't change with Ts, and that A and B are functions of Ts, and that B requires multiplication by G expressed in old units (i.e. per Ts0). Then our new T expressed in old units as a rate (i.e. dollars per Ts0) is just T = T(Ts)/Ts = ((1-A)/Ts)*H[n-1] + ((Ts-B)/Ts)*G[n], thus matching terms  we have  CTalt = (1-A)/Ts and DTalt = (Ts-B)/Ts, and because C/D = alpha2 (where C is Cm, not SIM's consumption output C), then alpha2alt = CTalt/DTalt. I think the idea came from SIM3. Also see SIM5. Talt is calculated in terms of CTalt and DTalt in the columns you need to scroll to the right to see. The alpha2alt isn't actually used for anything (I don't think). It's a given that this experiment will "work" and that Jason's Γ = 1 for all Ts (the purpose of the "experiment"), given that we could just directly calculate Talt from the other columns: namely Talt = G - ΔH/Ts where the Talt and G columns are expressed as dollars per Ts0. Thus ΔH/Ts is just ΔH expressed in those same units. But ultimately the alpha2alt I calculate here doesn't lead anywhere. In contrast, SIM6 is a successful approach: it adjusts both alpha1 and alpha2 (but not theta) to change sample periods without affecting steady state values (expressed in dollars per Ts0 and in dollars for the case of H) or the time constant Tc (expressed in Ts0). An update there of the date of this update explains why you adjust alpha1 and alpha2 but not theta to accomplish that. Given that Ts is changed by changing these parameters, then it's a given that Jason's Γ = 1 for all Ts in SIM6 as well. That's not really the case here in SIM2 because the resultant alpha2alt is not used (because there's probably not a way to make it useful). Thus really Jason's Γ 1 for all Ts here (i.e. using the method here to change sample  periods). And no, alpha2alt here is not equal to the alpha2 from SIM6. However, that's not to say that the SIM3 idea didn't ultimately lead to a greater success: see SIM9  and SIM10  for example.

The problem with my previous sheet was that the H curve changed when I changed the sample period (Ts). This was in spite of modifying the alpha2 parameter as per Ramanan's instructions.

I've changed tactics here: instead of scaling alpha2 (from a user provided alpha2_R parameter), I instead just use the User Parameters (in the green cells) directly without doing any modifications, (e.g. modifications such as scaling alpha2 by the ratio of change in the sample period (Ts/Ts_orig)).

What I do instead is record the original Ts setting (Ts_orig = 1 years). Then from the alpha1, alpha2 and theta parameters I derive continuous time parameters for the system. In particular I do the following:

Starting with the original discrete time system with A_orig and B_orig based on alpha1, alpha2 and theta, we have:

H[n+1] = A_orig*H[n] + B_orig*G[n+1]

From this I calculate the parameters for an equivalent continuous time system assuming G is held constant over each sample period:

dH(t)/dt = H'(t) = a_cont*H(t) + b_cont*G(t)

From this continuous time system, I re-descretize it with a new sample period (the user input sample period Ts in the green cells in the upper left). This gives me a new A and a new B, which I simply call A and B like in the old spreadsheet, and use them as follows in the results table:

H[n+1] = A*H[n] + B*G[n+1]

Here's how I do the calculations. First I assume the A_orig and B_orig resulted from a discetization of a continuous time system with an original sample period Ts = 1 year (Ts_orig = 1). Thus:

A_orig = expm(a_cont*Ts_orig) = exp(a_cont*Ts_orig)  (matrix exponentiation expm() reduces to scalar exponentiation exp() in this case, since the system has a single state).

B_orig = (A_orig - 1)*b_cont/a_cont  (to see why this is, look here).

I solve for a_cont and b_cont from the above, and then use these same formulas with a new sample period (Ts) to calculate A and B.

Some of these calculations are performed in the cells to the right of the chart, and the rest are in the cells for the A and B matrices respectively.

Matrices C and D (used for calculating the measurement vector [Y T YD C]' at sample n+1 by scaling H[n] and G[n+1] respectively), were not changed.

Now (you may notice) instead of H changing with Ts, ΔH changes with Ts. Is that a problem? I don't think so, because it SHOULD change with Ts: each Ts period of time now adds less to the total stock of money in the economy. H should NOT be changing with Ts, because the sample period should not affect the total stock of money accumulated in the economy at a fixed time in years.

Further information about calculating A and B as a function of Ts/Ts_orig: 
 DoΔ Δ
Starting with A_orig and B_orig and doing a little algebra reveals:

(1) A_orig = A(Ts_orig) = 1 - theta * alpha2 / (1 - alpha1 * (1 - theta))

(2) B_orig = B(Ts_orig) = 1 - theta / (1 - alpha1 * (1 - theta))

I write A(Ts_orig) and B(Ts_orig) because they are implicit functions of Ts_orig! Equations (1) and (2) were our starting point for A and B on the previous spreadsheet. Now what happens when we change Ts_orig to Ts?:

(3) A(Ts) = A_orig^(Ts/Ts_orig)

(4) B(Ts) = (B_orig/(A_orig - 1)) * (A(Ts) - 1) = -[(1 - alpha1*(1 - theta) - theta)/(theta*alpha2)] * (A(Ts) - 1)

Note: I didn't bother with that algebra in the spreadsheet above: I just used exp() and ln() according the the development above equation (1).


UPDATE 1: March 6, 2016 1:02 AM
I made a few changes to the plot to eliminate apparent sensitivity to Ts. Basically two things: I start the plot at t=0 years, but I start the table at -Ts. Secondly, instead of plotting H (which I leave in the main table), I plot H[n-1] (H delayed by one). H[n-1] takes the place of x[n] in a normal (canonical) discrete time state space formulation, so I suspected that I was plotting the wrong sequence against the others (off by one), and sure enough, once I changed that it helped quite a bit. The plot now shows near complete Ts invariance, as it should.

UPDATE 2: March 7, 2016 9:06 PM
Note that I resolved my problem with the old method of changing sample periods suggested by Ramanan on that older post, and I've updated it accordingly (read Update #2 there at the top). However that doesn't mean that the alpha2 scaling method is the same as this method: it is for small values of Ts, but try Ts >= 10 on this sheet and that one to see how they still differ.

UPDATE 3: March 8, 2016 12:33 AM
Note the  Γ below is Jason Smith's Γ from his March 2016 posts on SIM.


UPDATE 4: March 9, 2016, 5:21 PM
I made a new version called SIM3 posted here which combines the best of the version above (SIM2) with the old version (SIM), for the price of restricting the functional form of G.

UPDATE 5: March 13, 2016, 2:50 PM
I added some "alternative" parameters and variables to the spreadsheet. They don't affect any of the previous parameters or the plot, except that I've plotted T-alt to the plot. It differs slightly from T when Ts is changed. T-alt satisfies G&L's ΔH = G - T (actually ΔH = Ts*(G - T) since G and T are always rates expressed in dollars/(Ts_orig = Ts1 = 1 period = 1 year) in SIM2). Eventually there will be alternate versions of theta (thetaalt) and alpha1 (alpha1alt) as well. The idea is to eventually satisfy ALL of G&L's accounting equations, while maintaining sample period invariance, and matching G&L at Ts=1 (at least with G functions restricted to the form of their example). I have more information about it under the heading "SIM5" on a new post called SIM5.

UPDATE 6: March 17, 2016, 10:17 AM
As I explain in a new post (SIM6), and in a comment on Jason's blog, think Jason's Γ (above) represents both an approximate change to a new sample time (without actually changing the sample time values in the table), thus a time constant change without a change in steady states (what I explain in the comment). And, as I explain in the text of SIM6 below the spreadsheet) Jason's Γ as calculated above represents the difference between a sampled rate variable (like I treat G, Y, T, YD and C in the spreadsheet in this post) and an average rate calculated across a sample period interval (how they are treated in SIM and SIM6). As Ts gets larger, that difference grows.


  1. There are standard time units in economics. They are called 'days'. Every day starts, we do something and it ends. In the banking cycle intraday is a different process from inter-day and it is the same within business. Many of the fundamental processes are day based.

    Even the FX system has a cut off point of 10pm GMT where everybody has to settle.

    The DMO's operations are dynamic throughout the day and there are fixed points when gilts, etc issued and notified well in advance. Loads of co-ordination techniques used to reduce 'risk' and 'cost' over time in what is really an elaborate game.

    There's lots of good stuff on the DMO (UK Debt Management Office) site where it details that the DMO targets a weekly cash balance in their BoE account to allow for daily buffering:

    DMO site link

    So there is no continuous time. There are intra-day chunks and inter-days chunks with closing boundaries.

    As I said on Jason's blog:

    "In the UK BoE offers *unlimited intra-day overdrafts* and nobody operating those accounts needs to co-ordinate the balances until the end of the day. That's how it works.

    The Treasury spends asynchronously and with no direct relationship with the actions of the debt management office or the bank of england.

    What the debt management office does is maintain the reserve balance in the Treasury accounts at a set value - entirely for the sake of appearance.

    There is no control function here, it is pure politics. Government cheques don't bounce because there is nobody with the capacity to do that. If Treasury decides to spend, then the cheque clears.

    Dynamically it is better for all instrument sellers in the system for the overdrafts to be at their highest points *before* they start offering their own instruments into the market. Since that would guarantee the best bid.

    What I'm saying is that the DMO does what it can to get the best deal for the Treasury, and the best is likely to be achieved if it offers after it believes those overdrafts have expanded - since there will be more bids in the market and they would get a better price.

    I would expect serious treasury department in any organisation to know when they get the best price for their placements.

    So I would suggest it irrational behaviour for a treasury department of any organisation to pre-fund if they can get away with it."

    1. There may be no continuous time, but still the proper way to change sampling periods is:

      A2 = A1^(T2/T1)
      B2 = B1∙(A2-1)/(A1-1)

      Where, in this case A1 and B1 happen to be specified by:

      A1 = 1 - θ∙α2/(1 - α1 + θ∙α1)
      B1 = 1 - θ/(1 - α1 + θ∙α1)

      Scaling α2 by T2/T1 won't do it.

    2. Correction: scaling α2 by T2/T1 can be an acceptable approximation (see here), but the right way to do it (IMO) is still what I wrote in my above comment, as you can see by trying Ts=10 or larger on the spreadsheet above and at the old one at that link. They both work about the same for Ts < 3 or so, but large Ts values actually change the dynamics of the system using the α2 scaling method.

    3. Jason Smith makes an important observation about restrictions on alpha1 and alpha2 from G&L's book here.

  2. Looks like jason isn't letting my comment go through at his blog...

    It seems like you don't understand accounting units, For any stock or balance sheet items the unit is currency. For any flow or income statement the unit is currency/time period. Flows connect stocks at specific points in time.

    You should set up the Godley Lavoie flow and stock matrix's in your model, should make it more clear how they are connected.

    1. A H, thanks for you comment here. Did you see my response?

      In short I'll make what I expect to be minor changes to change the units of the exogenous input and the outputs so that inputs and outputs mimic an SFC model and appear "stock-flow consistent." There will be minor differences between the two (state space and SFC) because the compounding in state space models is truly sample period independent. For models exhibiting harmonic oscillation or exponential growth the differences will be more severe (I expect).

    2. ... that is "stock-flow consistent" over any particular sample period. Of course to compare them on a chart (say to check for sample time invariance) I sill think it makes sense to display them in fixed units which do not change with changing sample periods, as I already have it above. I expect those changes to be essentially cosmetic, and in no way to alter the system dynamics. And I already do have one SFC-style curve I plot: the ΔH curve (I explain in the post why it makes sense that it changes with changes in Ts, while the others do not).

    3. A H, OK, I took your suggestion, but with the old spreadsheet. Scaling alpha2 by T2/T1 works OK when I do that but does produce different results than the method above when T2 is too large. When T2 (the Ts parameter on both sheets) is close to 1 or less they are almost indistinguishable, but try setting Ts = 10 or 20 on both that sheet and this one and you can see that doing so on the old sheet changes the dynamics of the system profoundly (it starts to ring), while the one here continues to behave with the same dynamics (a damped rising exponential).

      Also, I think Jason did let (some?) of your comments through.

    4. Jason Smith makes an important observation about restrictions on alpha1 and alpha2 from G&L's book here.