Mata block nested in -capture- block, nested in loops

[David Mant]

Hi All,

I have a problem when a mata block is called within a -capture- block, which is itself nested within three Stata -for- loops.

The following simple code is to:

- iterate over various time series variables
- for each variable, estimate ARMA models for different lags
- store AIC and BIC for each model in a matrix
- show to the user the models that minimise the IC

mata: mata clear
local ar_max = 3
local ma_max = 3
// list of time series variables for analysis:
local varlist = "lond_r_rea east_r_rea seas_r_rea swes_r_rea"

foreach L_auth of local varlist {
    // Initialise empty matrix to store AIC, BIC of ARMA models
    mata {
        cnt = 0
        IC_`L_auth' = J((`ar_max')*(`ma_max'), 4, .)
    }

    // Loops over the ARMA models:
    forvalues i = 1(1)`ar_max' {
        forvalues j = 1(1)`ma_max' {
            // -capture- handles convergence issues in -arima- estimation
            capture noisily {
                display "Debug 0) This message shows the Local Authority: `L_auth'"
                arima `L_auth', ar(1/`i') ma(1/`j')
                estimates stats
                display "Debug A) This message should NOT be displayed if -arima- didn't converge: `L_auth'"
                // Saving AIC BIC (with ARMA order) in a single matrix, for each local authority:
                mata {
                    cnt = cnt + 1
                    IC_`L_auth'[cnt,] = (`i',`j',st_matrix("r(S)")[1,5..6])
                    IC_`L_auth'
                }
                display "Debug B) This message should NOT be displayed if -arima- didn't converge: `L_auth'"
            }
            display "Debug C) This message should display the Local Authority: `L_auth'"
        }
    }
    display "Debug D) This message should display the Local Authority: `L_auth'"
    
    // (the following part of the code is not problematic)
    mata {
        // Retrieve row where AIC = min and return the corresponding ARMA(p,q) order to Stata
        AIC_row = select((1::rows(IC_`L_auth')), (IC_`L_auth'[,3] :== colmin(IC_`L_auth'[,3])))
        st_numscalar("`L_auth'_AIC_ar", IC_`L_auth'[AIC_row,1])
        st_numscalar("`L_auth'_AIC_ma", IC_`L_auth'[AIC_row,2])

        // Retrieve row where BIC = min and return the corresponding ARMA(p,q) order to Stata
        BIC_row = select((1::rows(IC_`L_auth')), (IC_`L_auth'[,4] :== colmin(IC_`L_auth'[,4])))
        st_numscalar("`L_auth'_BIC_ar", IC_`L_auth'[BIC_row,1])
        st_numscalar("`L_auth'_BIC_ma", IC_`L_auth'[BIC_row,2])
    }
}
scalar list

The code works well, but I'm having issues in handling the loop when the -arima- estimation does not converge.
To handle such a scenario, I'm relying on the block:

capture noisily {
     ...
     mata {
     ...
     }
}

However, -capture- does not behave as I'd expect only when it contains a mata block.
In particular:

1 - if the -arima- estimation converges, no problem.

2 - if the -arima- estimation ends with a convergence error, my expectation is that the code that follows -arima- and is within the =capture- block is not executed. However, I'm finding two issues:
2a) the Stata code within the -capture- block that follows the mata block is executed anyway
2b) the local macro `L_auth' becomes blank

Below is the output under scenario 2. Note that:
- the message "Debug B)..." is displayed, even if it shouldn't
- the message "Debug D) ..." does not show the name of the local authority, even if it should.
(- the last error 3499 is a consequence of the previous one: the macro `L_auth' became unexpectedly blank so it cannot find the matrix)

...
Iteration 247: log likelihood = 860.05825 (backed up)
flat log likelihood encountered, cannot find uphill direction
Debug B) This message should NOT be displayed if -arima- didn't converge: lond_r_rea
Debug C) This message should display the Local Authority: lond_r_rea
Debug D) This message should display the Local Authority:
<istmt>: 3499 IC_ not found
r(3499);

end of do-file

r(3499);

.

And below is the output still after an -arima- convergence error, but I've now simply deleted the mata block contained within the ARMA loops.
Note that the output is now perfectly sensible to me:
- the message "Debug B)..." is not displayed. Correct.
- the message "Debug D) ..." shows the name of the local authority. Correct.
(- the error <istmt> at the very end is expected as I've deleted the aforementioned mata block, so no problem)

Iteration 247: log likelihood = 860.05825 (backed up)
flat log likelihood encountered, cannot find uphill direction
Debug C) This message should display the Local Authority: lond_r_rea
Debug D) This message should display the Local Authority: lond_r_rea
st_numscalar(): 3204 matrix found where scalar required
<istmt>: - function returned error
r(3204);

end of do-file

r(3204);

.

In a nutshell, it is as if the mata block causes issues within the loops / blocks.
Any suggestion would be greately appreciated.
Thank you.

[Niels Henrik Bruun]

Hi
I would suggest that You separate Stata and Mata code.
Build Mata funktions to call from Stata.
Consider Mata loops instead of Stata ditto.

Kind regards
nhb

[David Mant]

Hi Niels,

thank you for your swift response.
I was hoping there was a way to keep the mata block as is, i.e. using {...}, rather than having to build functions.
The reason is to improve readability (I have other similar instances in other parts of the code).

I have now stripped the code above down to the bone, to isolate the issue:

mata: mata clear
local varlist = "lond_r_rea east_r_rea seas_r_rea swes_r_rea"

foreach L_auth of local varlist {
    forvalues i = 1(1)3 {
        forvalues j = 1(1)3 {
            capture noisily {
                display "Debug 0) This message shows the Local Authority: `L_auth'"
                arima `L_auth', ar(1/`i') ma(1/`j')
                display "Debug A) This message should NOT be displayed if -arima- didn't converge: `L_auth'"
                mata {
                    a = 10
                    b = 20
                 }
                display "Debug B) This message should NOT be displayed if -arima- didn't converge: `L_auth'"
            }
            display "Debug C) This message should display the Local Authority: `L_auth'"
        }
    }
    display "Debug D) This message should display the Local Authority: `L_auth'"
}

If I run the code when -arima- does not converge, I get this output:

...
Iteration 247: log likelihood =  860.05825  (backed up)
flat log likelihood encountered, cannot find uphill direction
Debug B) This message should NOT be displayed if -arima- didn't converge: lond_r_rea
Debug C) This message should display the Local Authority: lond_r_rea
Debug D) This message should display the Local Authority:
} is not a valid command name
r(199);

end of do-file

r(199);

I also flag that the mata block is executed, even if it should not.
It is as if Stata is not able to correctly interpret the nested curly brackets of my code above.
More precisely, it is as if Stata added a brace right before the mata block.
In fact, the output above would be perfectly consistent with adding a curly bracket where shown below:

mata: mata clear
local varlist = "lond_r_rea east_r_rea seas_r_rea swes_r_rea"

foreach L_auth of local varlist {
    forvalues i = 1(1)3 {
        forvalues j = 1(1)3 {
            capture noisily {
                display "Debug 0) This message shows the Local Authority: `L_auth'"
                arima `L_auth', ar(1/`i') ma(1/`j')
                display "Debug A) This message should NOT be displayed if -arima- didn't converge: `L_auth'"

            } // IT IS AS IF STATA ADDED THIS BRACE HERE BY ITSELF

            mata {
                A = 10
                A
                B = 20
                B
            }
            display "Debug B) This message should NOT be displayed if -arima- didn't converge: `L_auth'"
        }
        display "Debug C) This message should display the Local Authority: `L_auth'"
    }
}
display "Debug D) This message should display the Local Authority: `L_auth'"
}

That "automatically added" brace would explain:

- why mata code is executed
- why the message "Debug B) ..." is shown
- why the message "Debug D) ..." does not show "lond_r_rea" (i.e. the current value of `L_auth' in the -foreach- loop)
- why this error is ultimately returned: "} is not a valid command name", as now the last brace would be redundant.

Am I doing some mistakes in calling the mata block, or is this really a Stata parsing/compiling issue?
Many thanks.

[Niels Henrik Bruun]

Hi David
I'm sorry, that I misread your post.
If I try to reproduce your case with Mata code following a Stata error, things works as expected, i.e. the Mata isn't run:

. mata a = 0; b = 0

. forvalues r = -1/1 {
  2.         capture noisily {
  3.                 display "nhb `r'"
  4.                 mata {
  5.                         a++
  6.                         "a", strofreal( (a,b) )
  7.                 }
  8.                 if `r' == 0 error 134
  9.                 mata {
 10.                         b++
 11.                         "b", strofreal( (a,b) )
 12.                 }
 13.         }
 14. }
nhb -1
       1   2   3
    +-------------+
  1 |  a   1   0  |
    +-------------+
       1   2   3
    +-------------+
  1 |  b   1   1  |
    +-------------+
nhb 0
       1   2   3
    +-------------+
  1 |  a   2   1  |
    +-------------+
too many values
nhb 1
       1   2   3
    +-------------+
  1 |  a   3   1  |
    +-------------+
       1   2   3
    +-------------+
  1 |  b   3   2  |
    +-------------+

[Niels Henrik Bruun]

However, if your Mata block is outside your capture block, it should be run. Maybe you should set your capture blocks differently.

[David Mant]

Thank you Niels. I have tried many (small) coding variations, but still the nested mata block was being executed somehow...
I ended up doing things in a completely different way to bypass the issue.
All the best

[David Mant]

I have come across the aforementioned issue in a more general (and simpler to discuss) form.
I thought it would have been more practical to start a new thread, as here I've depicted the problem in a more complex/specific way than it actually is (as I've now realised).
Thank you.