Statalist The Stata Listserver


[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]

Re: st: mata: why this loop is so slow


From   wgould@stata.com (William Gould, Stata)
To   statalist@hsphsun2.harvard.edu
Subject   Re: st: mata: why this loop is so slow
Date   Tue, 30 May 2006 12:41:43 -0500

AbdelRahmen <uaquap@gnet.tn> reported

> I've written this portion of code (see below) and it works
> perfectly the problem is the execution time which is about 15
> mn. Under SAS (IML) my co-author code the same thing and the
> execution time is 1 (one) minute on 35900 observations. is
> there a problem with my programmation Style?

I replied, AbdelRahmen replied, and I now have run his code on my 
computer.  

It ran is 0.84 seconds on my computer.

Perhaps I have done something wrong.  I concocted phony values, but 
everything I think is as AbdelRahmen describes it.  I wonder how 
fast the attached do-file runs on his computer.

-- Bill
wgould@stata.com

File try.do follows:
------------------------------------------------------------------------------
cscript

set memory 30m
set rmsg on
local OBS 35900

set obs `OBS'
set seed 1234

forvalues i=1(1)46 {
        gen v`i' = uniform()
}
gen vprobit = uniform()
gen dmills  = uniform()
gen selection = uniform()
gen t = mod(_n-1, 5) + 1

local instrume v1  v2  v3  v4  v5  v6  v7  v8  v9  v10 ///
               v11 v12 v13 v14 v15 v16 v17 v18 v19 v20 ///
               v21 v22 v23 v24 v25 v26 v27 v28 v29 v30 ///
               v31 v32 v33 v34 v35 v36 v37 v38 v39 v40 ///
               v41 v42 v43 v44 v45 v46 

local vprobit "vprobit"

mata:

w = st_data(., tokens(st_local("instrume")))    // matrix
q = st_data(., tokens(st_local("vprobit")))     // matrix
dmills=st_data(., "dmills")                     // colvector
s=st_data(., "selection")                       // colvector
t=st_data(., "t")                               // colvector

chw=cols(w)                                     // scalar
cq=cols(q)                                      // scalar

h = uniform(`OBS', 1)

n = nt = `OBS'
T  = 21
beta = uniform(46,1)

F = J(1,21,0)
for (i=1; i<=nt; i++) {
        jacobi=J(chw, (cq*T) ,0)
        jacobi[(chw - T + t[i]), (cq*t[i] -cq+1)..cq*t[i]] = 
                                                   q[i,.]:*dmills[i,.]
        F = F + (s[i,.]:*h[i,.]')*beta'*jacobi
}

F=(1/n)*F

end
------------------------------------------------------------------------------
<end>
*
*   For searches and help try:
*   http://www.stata.com/support/faqs/res/findit.html
*   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   |   What's new   |   Site index