/*************************************

mataclass.do

written 10/2007 by Russell Dimond

Examples from
http://ssc.wisc.edu/sscc/pubs/4-26.htm

see also ex1.do and ex2.do

************************************/

mata
mata clear
end
clear
log using mataclass.log, replace

/* Basics */
mata
2+2
x=2+2
x

/* Column and Row Join Operators */
1,2
1\2
1,2\3,4

/* Errors */
1,2\3\4
1,\3\4

/* Parentheses */
1,2\3,4
(1,2\3,4)
(1,2)\(3,4)

/* Range Operators */
1..3
1::3
1.2..3.5

/* Variables */
x=3,4
y=5,6
z=(1,2)\x\y
z
/* Arithmatic Operators */
x=(1,2)\(3,4)
y=(1,2)\(3,4)
x*y
x:*y

/* Logical Operators */
x==y
x>y
x<y
z=(1,2)\(3,5)
x==z
x>z
x<z
x<=z
x:==y
x:==z
x:<z
x:==2 
x:==(1,2)
x:==(1\3)

/* Transpose Operator */
x'

/*Subscripts */
z=x[1,1]
z
x[1,1]=5
x
y=1,2,3
y[2]
y=1\2\3
y[2]
x[.,1]
x[1,]
a=(1,2,3)\(4,5,6)\(7,8,9)
a[(1,3),.]
a[(3,2,1),(3,2,1)]
a[(3,1,2,2,1,3),.]
b=(3,1,1)
a[.,b]

/* Subscripting Ranges */
x1=1..10
x2=(1::10)*10
x=x1[(1,1,1,1,1,1,1,1,1,1),.]+x2[.,(1,1,1,1,1,1,1,1,1,1)]
x
x[(4..7),(3..8)]
x[(4,5,6,7),(3,4,5,6,7,8)]
x[(4..7),(3..8)]
x[|4,3\7,8|]
x[|.,.\.,.|]
x[|3,3\3,.|]

/* Getting your Data from Stata to Mata */
end
sysuse auto
l in 1 
mata

/* st_data */
_st_data(1,2)
_st_sdata(1,1)
_st_data(1,st_varindex("mpg"))
st_data(.,.)
st_data((1,3),2)
st_data(1,(2,4))
st_data(1,("price","rep78"))
st_sdata(.,1)
x=st_data(1,(2,4))
x
st_sdata(.,1,"foreign")

/* st_view */
st_view(x,.,.)
x
st_view(n=0,1,(2,4))
n

/* Getting your Data from Mata to Stata */
st_view(mpg=0,.,"mpg")
st_view(weight=0,.,"weight")
pmg=weight:*mpg
st_addvar("long","pmg")
st_store(.,"pmg",pmg)
end
by foreign: sum pmg

/* Hierarchical Data */
use hh,replace
l
mata
mata matuse hh
hh
ind
hh[ind[2,3],1] 
x=ind,hh[ind[.,3],.]

/* Matrix Functions */
I(3)
I(4,3)
J(3,3,0)
J(2,3,"a")
e(1,3)
e(1,3)\e(2,3)\e(3,3)
uniform(5,5)
x=(1,2)\(3,4)
y=makesymmetric(x)
y
_makesymmetric(x)
x

/* Sorting */
x=(2,1)\(1,3)\(1,2)
sort(x,1)
sort(x,(1,2))
x
_sort(x,(1,2))
x

/* Size of Matrices */
rows(x)
cols(x)
length(x)
y=(1,2,.)\(4,.,6)\(7,8,.)
missing(y)
rowmissing (y)
colmissing(y) 

/* Descriptive Statistics */
sum(x)
rowsum(x)
colsum(x)
max(x)
min(x)
rowmax(x)
colmin(x)
minmax(x)
colminmax(x)
mean(x)
rowsum(x)/cols(x)
variance(x)
correlation(x)

/* Matrix Characteristics */
x=(1,2)\(3,4)
diagonal(x) 
trace(x)
det(x)
rank(x)

/* Solvers, Inverters and Decomposition */
a=uniform(5,5)
a=a'*a
cholesky(a) 
cholesky(a)*cholesky(a)'
b=uniform(5,1)
cholsolve(a,b)
a*cholsolve(a,b)
cholinv(a)
a*cholinv(a)
a=uniform(5,5)
lud(a,L=0,U=0,p=0)
(L*U)[p,.]
a*lusolve(a,b)
a*luinv(a)

/* Linear Regression */
end
sysuse auto
gen touse=(price!=. & mpg!=. & rep78!=. & weight!=. & foreign!=.)
mata
st_view(x=0,.,("mpg","rep78","weight","foreign"), "touse")
st_view(y=0,.,("price"), "touse")
x=x,J(rows(x),1,1)
b=invsym(x'*x)*x'*y
e=y-x*b
v=(e'*e)/(rows(x)-cols(x))*invsym(x'*x)
se=sqrt(diagonal(v))
t=b:/se
p=2*ttail(rows(x)-cols(x),abs(t))
b,se,t,p
end
reg price mpg rep78 weight foreign

/* Loops */
mata
x=0
while (x<=5)
{
	x
	x++
}
while (x<=5) x++
x

do
{	
	x++
} while (x<=5)
x

for (x=0; x<=5; x++)
{
	x
} 

for(; x<=10;) {
	x
	x++
}

m=J(5,5,0)
for(i=1; i<=rows(m); i++) {
	for(j=1; j<=cols(m); j++) {
		m[i,j]=i+(j-1)*cols(m)
	}
} 
m

/* Writing Your Own Functions */
function doubleAndSum(x) {
x=x*2
return(sum(x))
}

m=I(3)
doubleAndSum(m)
m
doubleAndSum(I(3))
doubleAndSum(1..3)



/* Pointers */
mata
x=(1..3)\(4..6)\(7..9)
p=&x
*p
(*p)[2,2]
p=&J(5,5,0)

x=J(5,5,NULL)
for (i=1; i<=rows(x); i++)
{
	for (j=1; j<=cols(x); j++)
	{
		x[i,j]=&J(5,1,0)
	}
}

function put(val,x,i,j,k)
{
	/* Usage: value to put, matrix to put it in, i, j, k to put it at. */
	(*(x[i,j]))[k,1]=val
}

function get(x,i,j,k)
{
	/* Usage: matrix to get from, i, j, k of value to get. */ 
	return((*(x[i,j]))[k,1])
}

put(3,x,5,5,5)
get(x,5,5,5) 

p=&doubleAndSum()
(*p)(I(3))

/* Mata's Optimizer */
function f(todo,x,y,g,H)
{
	y=-x^4
} 

function g(todo,x,y,g,H)
{
	y=-x^4
	if (todo>=1)
	{
		g=-4*x^3
		if (todo==2)
		{
			H=-12*x^2
		}
	}
}

function h(todo,x,y,g,H)
{
	y=-x[1]^4-x[2]^4
}

function i(todo,x,y,g,H)
{
	y=J(2,1,.) 
	y[1]=-x[1]^4
	y[2]=-x[2]^4 
} 

s=optimize_init()
optimize_init_evaluator(s,&f())
optimize_init_params(s,-100)
optimize(s)
xmax=optimize(s)

s=optimize_init()
optimize_init_evaluator(s,&g())
optimize_init_evaluatortype(s,"d2")
optimize_init_params(s,-100)
optimize(s)

s=optimize_init()
optimize_init_evaluator(s,&h())
optimize_init_params(s,(-100,100))
optimize(s)

s=optimize_init()
optimize_init_evaluator(s,&i())
optimize_init_evaluatortype(s,"v0")
optimize_init_params(s,(-100,100))
optimize(s)

end
log close