Programming an estimation command in Stata: Mata 101
I introduce Mata, the matrix programming language that is part of Stata.
This is the eleventh post in the series Programming an estimation command in Stata. I recommend that you start at the beginning. See Programming an estimation command in Stata: A map to posted entries for a map to all the posts in this series.
Meeting Mata
Mata is a matrix programming language that is part of Stata. Mata code is fast because it is compiled to object code that runs on a virtual machine; type help m1_how for details.
The easiest way to learn Mata is to use it. I begin with an interactive session. (You might find it useful to type along.)
Example 1: A first interactive Mata session
. mata:
------------------------------------------------- mata (type end to exit) ------
: X = J(3, 4, 5)
: X
1 2 3 4
+-----------------+
1 | 5 5 5 5 |
2 | 5 5 5 5 |
3 | 5 5 5 5 |
+-----------------+
: w = (1::4)
: w
1
+-----+
1 | 1 |
2 | 2 |
3 | 3 |
4 | 4 |
+-----+
: v = X*w
: v
1
+------+
1 | 50 |
2 | 50 |
3 | 50 |
+------+
: v'
1 2 3
+----------------+
1 | 50 50 50 |
+----------------+
: end
--------------------------------------------------------------------------------
Typing mata: causes Stata to drop down to a Mata session. Typing end ends the Mata session, thereby popping back up to Stata. The dot prompt . is Stata asking for something to do. After you type mata:, the colon prompt : is the Mata compiler asking for something to do.
Typing X = J(3, 4, 5) at the colon prompt causes Mata to compile and execute this code. J(r, c, v) is the Mata function that creates an r\(\times\)c matrix, each of whose elements is v. The expression on the right-hand side of the assignment operator = is assigned to the symbol on the left-hand side.
Typing X by itself causes Mata to display what X contains, which is a 3\(\times\)4 matrix of 5s. Unassigned expressions display their results. Type help m2_exp for details about expressions.
Typing w = (1::4) causes Mata to use the column range operator to create the 4\(\times\)1 column vector that was assigned to w and displayed when I typed w by itself. Type help m2_op_range for details and a discussion of the row range operator.
Typing v = X*w causes Mata to assign the matrix product of X times w to v, which I subsequently display. I then illustrate that ‘ is the transpose operator. Type help m2_exp, marker(remarks7) for a list of operators.
Again, typing end ends the Mata session.
In almost all the work I do, I extract submatrices from a matrix.
Example 2: Extracting submatrices from a matrix
. mata:
------------------------------------------------- mata (type end to exit) ------
: rseed(1234)
: W = runiform(4,4)
: W
1 2 3 4
+---------------------------------------------------------+
1 | .9472316166 .0522233748 .9743182755 .9457483679 |
2 | .1856478315 .9487333737 .8825376215 .9440776079 |
3 | .0894258515 .7505444902 .9484983174 .1121626508 |
4 | .4809064012 .9763447517 .1254975307 .7655025515 |
+---------------------------------------------------------+
: v = (2, 4)
: u = (1\ 3)
: v
1 2
+---------+
1 | 2 4 |
+---------+
: u
1
+-----+
1 | 1 |
2 | 3 |
+-----+
: W[u, v]
1 2
+-----------------------------+
1 | .0522233748 .9457483679 |
2 | .7505444902 .1121626508 |
+-----------------------------+
: W[| 1,1 \ 3,3 |]
1 2 3
+-------------------------------------------+
1 | .9472316166 .0522233748 .9743182755 |
2 | .1856478315 .9487333737 .8825376215 |
3 | .0894258515 .7505444902 .9484983174 |
+-------------------------------------------+
: end
--------------------------------------------------------------------------------
I use rseed() to set the seed for the random-number generator and then use runiform(r,c) to create a 4\(\times\)4 matrix uniform deviates, which I subsequently display.
Next, I use the row-join operator , to create the row vector v and I use the column-join operator \ to create the column vector u. Type help m2_op_join for details.
Typing W[u,v] extracts from W the rows specified in the vector u and the columns specified in the vector v.
I frequently extract rectangular blocks defined by a top-left element and a bottom-right element. I illustrate this syntax by typing
W[| 1,1 \ 3,3 |]
In detail, [| opens a range-subscript extraction, 1,1 is the address of the top-left element, \ separates the top-left element from the bottom-right element, 3,3 is the address of the bottom-right element, and |] closes a range-subscript extraction. Type help m2_subscripts for details.
Ironically, when I am doing matrix programming, I frequently want the element-by-element operator instead of the matrix operator. Preface any matrix operator in Mata with a colon (:) to obtain the element-by-element equivalent.
Example 3: Element-wise operators
. mata:
------------------------------------------------- mata (type end to exit) ------
: W = W[| 2,1 \ 4,4 |]
: W
1 2 3 4
+---------------------------------------------------------+
1 | .1856478315 .9487333737 .8825376215 .9440776079 |
2 | .0894258515 .7505444902 .9484983174 .1121626508 |
3 | .4809064012 .9763447517 .1254975307 .7655025515 |
+---------------------------------------------------------+
: v = .1*(4::6)
: v
1
+------+
1 | .4 |
2 | .5 |
3 | .6 |
+------+
: v:*W
1 2 3 4
+---------------------------------------------------------+
1 | .0742591326 .3794933495 .3530150486 .3776310432 |
2 | .0447129257 .3752722451 .4742491587 .0560813254 |
3 | .2885438407 .585806851 .0752985184 .4593015309 |
+---------------------------------------------------------+
: v'*W
1 2 3 4
+---------------------------------------------------------+
1 | .4075158991 1.340572446 .9025627257 .8930138994 |
+---------------------------------------------------------+
: end
--------------------------------------------------------------------------------
I extract the bottom four rows of W, store this matrix in W, and display this new W. I then create a row-wise conformable vector v, perform element-wise multiplication of v across the columns of W, and display the result. I cannot type v*W because the 3\(\times\)1 v is not conformable with the 3\(\times\)3 W. But I can, and do, type v’*W because the 1\(\times\)3 v’ is conformable with the 3\(\times\)3 W.
Example 4 uses an element-wise logical operator.
Example 4: Element-wise logical operator
. mata:
------------------------------------------------- mata (type end to exit) ------
: W :< v
1 2 3 4
+-----------------+
1 | 1 0 0 0 |
2 | 1 0 0 1 |
3 | 1 0 1 0 |
+-----------------+
: end
--------------------------------------------------------------------------------
I display the result of comparing the element-wise conformable v with W. Type help m2_op_colon for details.
Stata data in Mata
The Mata function st_data() creates a Mata matrix containing a copy of the data from the Stata dataset in memory. The Mata function st_view() creates a Mata view of the data in the Stata dataset in memory. Views act like matrices, but there is a speed-space tradeoff. Copies are fast at the cost of using twice as much memory. Views are slower, but they use little extra memory.
Copying the data from Stata into Mata doubles the memory used, but the values are stored in Mata memory. Every time a Mata function asks for a value from a matrix, it finds it immediately. In contrast, a view of the data in Stata barely increases the memory used, but the values are in Stata memory. Every time a Mata function asks for a value from a view, it finds a sign telling it where in Stata to get the value.
Example 5: Data from Stata into Mata
. sysuse auto
(1978 Automobile Data)
. list mpg headroom trunk rep78 turn foreign in 1/3 , nolabel
+-------------------------------------------------+
| mpg headroom trunk rep78 turn foreign |
|-------------------------------------------------|
1. | 22 2.5 11 3 40 0 |
2. | 17 3.0 11 3 40 0 |
3. | 22 3.0 12 . 35 0 |
+-------------------------------------------------+
. mata:
------------------------------------------------- mata (type end to exit) ------
: Y = st_data(., "mpg headroom trunk")
: st_view(X=., ., "rep78 turn foreign")
: V = Y,X
: V[| 1,1 \ 3,6 |]
1 2 3 4 5 6
+-------------------------------------+
1 | 22 2.5 11 3 40 0 |
2 | 17 3 11 3 40 0 |
3 | 22 3 12 . 35 0 |
+-------------------------------------+
: X[3,1] = 7
: X[| 1,1 \ 3,3 |]
1 2 3
+----------------+
1 | 3 40 0 |
2 | 3 40 0 |
3 | 7 35 0 |
+----------------+
: end
--------------------------------------------------------------------------------
. list rep78 turn foreign in 1/3 , nolabel
+------------------------+
| rep78 turn foreign |
|------------------------|
1. | 3 40 0 |
2. | 3 40 0 |
3. | 7 35 0 |
+------------------------+
After I list out the first three observations on six variables in the auto dataset, I drop down to Mata, use st_data() to put a copy of all the observations on mpg, headroom, and trunk into the Mata matrix Y, and use st_view() to create the Mata view X on to all the observations on rep78, turn, and foreign.
After row-joining Y and X to create V, I display the first 3 rows of V. Note that the third observation on rep78 is missing and that Mata matrices and views can contain missing values.
Changing the value of an element in a view changes the data in Stata. I illustrate this point by replacing the (3,1) element of the view X with 7, displaying the first three rows of the view, and listing out the first three observations on rep78, turn, and foreign.
Copying matrices between Mata and Stata
The Mata function st_matrix() puts a copy of a Stata matrix into a Mata matrix, or it puts a copy of a Mata matrix into a Stata matrix. In example 6, V = st_matrix("B") puts a copy of the Stata matrix B into the Mata matrix V.
Example 6: Creating a copy of a Stata matrix in a Mata vector
. matrix B = (1, 2\ 3, 4)
. matrix list B
B[2,2]
c1 c2
r1 1 2
r2 3 4
. mata:
------------------------------------------------- mata (type end to exit) ------
: V = st_matrix("B")
: V
1 2
+---------+
1 | 1 2 |
2 | 3 4 |
+---------+
: end
--------------------------------------------------------------------------------
In example 7, st_matrix("Z", W) puts a copy of the Mata matrix W into the Stata matrix Z.
Example 7: Creating a copy of a Mata matrix in a Stata vector
. mata:
------------------------------------------------- mata (type end to exit) ------
: W = (4..6\7..9)
: W
1 2 3
+-------------+
1 | 4 5 6 |
2 | 7 8 9 |
+-------------+
: st_matrix("Z", W)
: end
--------------------------------------------------------------------------------
. matrix list Z
Z[2,3]
c1 c2 c3
r1 4 5 6
r2 7 8 9
Strings
Mata matrices can be string matrices.
In my work, I frequently have a list of variables in a string scalar that is easier to work with as a string vector.
Turning a string scalar list into a string vector
. mata:
------------------------------------------------- mata (type end to exit) ------
: s1 = "price mpg trunk"
: s1
price mpg trunk
: s2 = tokens(s1)
: s2
1 2 3
+-------------------------+
1 | price mpg trunk |
+-------------------------+
: end
--------------------------------------------------------------------------------
I use tokens() to create the string vector s2 from the string vector s1.
Flow of control
Mata has constructs for looping over a block of code enclosed between curly braces or only executing it if an expression is true.
I frequently use the for() construction to loop over a block of code.
mata:
for(i=1; i<=3; i=i+1) {
i
}
end
In this example, I set i to the initial value of 1. The loop will continue as long as i is less than or equal to 3. Each time through the loop, the block of code enclosed between the curly braces is executed, and 1 is added to the current value of i. The code block displays the value of i. Example 9 illustrates these points.
Example 9: A for loop
. mata:
------------------------------------------------- mata (type end to exit) ------
: for(i=1; i<=3; i=i+1) {
> i
> }
1
2
3
: end
--------------------------------------------------------------------------------
Sometimes, I want to execute a block of code as long as a condition is true, in which case I use a while loop, as in code block 2 and example 10.
i = 7
while (i>5) {
i
i = i - 1
}
I set i to 7 and repeat the block of code between the curly braces while i is greater than 5. The block of code displays the current value of i, then subtracts 1 from i.
Example 10: A while loop
. mata:
------------------------------------------------- mata (type end to exit) ------
: i = 7
: while (i>5) {
> i
> i = i - 1
> }
7
6
: end
--------------------------------------------------------------------------------
The if construct only executes a code block if an expression is true. I usually use the if-else construct that executes one code block if an expression is true and another code block if the expression is false.
. mata:
-------------------------------------------- mata (type end to exit) ---
: for(i=2; i<10; i=i+5) {
> i
> if (i<3) {
> "i is less than 3"
> }
> else {
> "i is not less than 3"
> }
> }
2
i is less than 3
7
i is not less than 3
: end
-------------------------------------------------------------------------
One-line calls to Mata
I frequently make one-line calls to Mata from Stata. A one-line call to Mata causes Stata to drop to Mata, compile and execute the line of Mata code, and pop back up to Stata.
Example 12: One-line calls to Mata
. mata: st_matrix("Q", I(3))
. matrix list Q
symmetric Q[3,3]
c1 c2 c3
r1 1
r2 0 1
r3 0 0 1
In example 12, I use the one-line call to Mata mata: st_matrix("Q", I(3)) to put a copy of the Mata matrix returned by the Mata expression I(3) into the Stata matrix Q. After the one-line call to Mata, I am back in Stata, so I use matrix list Q to show that the Stata matrix Q is a copy of the Mata matrix W.
Done and undone
I used an interactive session to introduce Mata, the matrix programming language that is part of Stata.
In the next post, I show how to define Mata functions.