Bookmark and Share

Notice: On March 31, it was announced that Statalist is moving from an email list to a forum. The old list will shut down on April 23, and its replacement, statalist.org is already up and running.


[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: st: looping in MATA


From   wgould@stata.com (William Gould, StataCorp LP)
To   statalist@hsphsun2.harvard.edu
Subject   Re: st: looping in MATA
Date   Fri, 30 Apr 2010 09:37:04 -0500

> I need to create a rolling loop in stata mata language, I am using
> some simple matrix algebra and it is easier to do in Mata than in
> matrix functions.
> 
> I need to do something like this (within a do file):
> *****************************************
> mata
> 
> foreach y in 1970/2000 {
> 
> x`y' = st_data(., (" g`y'", " g`y+1'", " g`y+2'"," g`y+3'", " g`y+4'"), 0)
> 
> r=rows(x`y')
> vc=x`y'*x`y''
> vc=vc/25
> var`y'=sqrt(trace(vc)/rows(vc))
> tr=vc-diag(vc)
> cov`y'=sum(tr)
> cov`y'=sqrt(cov`y'/(r^2-r)/2)
>
> .end
> *************************
>
> Ideally I would also like to save these cov`y' and var`y' data into a matrix.

Here is a solution:

        mata:

        result = J(0, 3, .)
        for (y=1970; y<=2000; y++) {
                v0 = sprintf("g%f", y)
                v1 = sprintf("g%f", y+1)
                v2 = sprintf("g%f", y+2)
                v3 = sprintf("g%f", y+3)
                v4 = sprintf("g%f", y+4)

                x  = st_data(., (v0, v1, v2, v3, v4))
                r  = rows(x)
                vc = (x*x')/25        // I think 25 should be r
                var= sqrt(trace(vc)/25)
                tr = vc - diag(vc)
                cov= sum(tr)
                cov= sqrt(cov/(r^2-r)/2)

                result = result \ (y,var,cov)
        }
        result

        end

N.B., I did not use macros!  As Nick Cox <n.j.cox@durham.ac.uk> Cox said in
his reply, "The key principle here is that local macro references are not part
of Mata. So, they must be interpreted before Mata sees them."  Nick went from
there to suggesting a a "clumsy" (his word) solution that called Mata line by
line.

I eliminated the macros altogether.  Whenever you're thinking macros and 
Mata simultaneously, you are thinking wrong, and the solution is not to 
try to figure out how to make the macros work, but instead to think about 
how to eliminate the macros altogether.

That was easy.

Olga thought, "I need to use use variables named g`y', g`=y+1', ..., where 
y = 1970, 1971, ...  That's the right way to think in Stata, but not 
in Mata.  Whereever you use macros in Stata, you need yet another variable
in Mata.  Variables in Mata are used for everything.  I created a sequence
of Mata string variables to hold the names I wanted to access:

                v0 = sprintf("g%f", y)
                v1 = sprintf("g%f", y+1)
                v2 = sprintf("g%f", y+2)
                v3 = sprintf("g%f", y+3)
                v4 = sprintf("g%f", y+4)

I used sprintf() to assemble the name, but I could just as well have coded

                v0 = "g" + strofreal(y) 
                v1 = "g" + strofreal(y+1) 
                v2 = "g" + strofreal(y+2) 
                v3 = "g" + strofreal(y+3) 
                v4 = "g" + strofreal(y+4) 

Code it however you wish, but understand, variable y is a real scalar, 
and I'm forming variable names and stroging them in variables v0, v1, ...,
v4, which are string scalars.  Once I had the variables containing the 
names of the Stata variables, I was ready to form matrix x:

                x  = st_data(., (v0, v1, v2, v3, v4))

By the way, I collected the results together in a single matrix, just 
as Olga requested.

By the way, I think Olga needs to check that she has the correct formulas 
for the calculation.  There is a 25 in the code that I believe should be 
r, and I note that the calculation as written is scale dependent.  To 
wit, I tested the above code on some literally (psuedo) random data.
Then I took the data and added 10 to to every variable.  Answers changed.

The scale dependencies can be traced to the line, 

                vc = (x*x')/25        // I think 25 should be r

I suspect Olga may be thinking of matrix x as a matrix of values recorded
in deviation from the mean form.  If so, changing absolute x into deviation 
form is easy enough: 

               x = x :- (colsum(x):/rows(x)) 

or, if you prefer:

               mean = colsum(x) :/ rows(x)
               x    = x :- mean

-- Bill
wgould@stata.com
*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/


© Copyright 1996–2014 StataCorp LP   |   Terms of use   |   Privacy   |   Contact us   |   Site index