-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathmvdom.ado
149 lines (90 loc) · 3.91 KB
/
mvdom.ado
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
*! mvdom version 1.2.1 12/20/2024 Joseph N. Luchman
version 16
program define mvdom, eclass
syntax varlist(min = 2) [if] [aw fw], dvs(varlist min=1) [epsilon pxy] //epsilon is a hidden option
tempname canonmat
tempvar touse
gettoken dv ivs: varlist
quietly generate byte `touse' = 1 `if'
quietly replace `touse' = 0 if missing(`touse')
if strlen("`epsilon'") {
mata: eps_ri_mv("`dv' `dvs'", "`ivs'", "`touse'")
}
else {
if strlen("`pxy'") {
quietly correlate `dv' `dvs' `ivs' [`weight'`exp'] if `touse'
local dvnum: word count `dv' `dvs'
matrix `canonmat' = r(C)
matrix `canonmat' = trace( ///
invsym(`canonmat'[1..`:word count `dv' `dvs'', 1..`:word count `dv' `dvs''])* ///
`canonmat'[`=`:word count `dv' `dvs''+1'..., 1..`:word count `dv' `dvs'']'* ///
invsym(`canonmat'[`=`:word count `dv' `dvs''+1'..., `=`:word count `dv' `dvs''+1'...])* ///
`canonmat'[`=`:word count `dv' `dvs''+1'..., 1..`:word count `dv' `dvs''] ///
)
}
else {
quietly _canon (`dv' `dvs') (`ivs') [`weight'`exp'] if `touse'
matrix `canonmat' = e(ccorr)
}
ereturn post, esample(`touse')
if strlen("`pxy'") ereturn scalar r2 = `canonmat'[1, 1]/`:word count `dv' `dvs''
else ereturn scalar r2 = `canonmat'[1, 1]^2
ereturn local title "Multivariate regression"
}
end
/*Mata function to execute epsilon-based relative importance with mvdom*/
version 16
mata:
mata set matastrict on
void eps_ri_mv(string scalar dvlist, string scalar ivlist, string scalar touse)
{
/*object declarations*/
real matrix X, Y, L, R, Lm, L2, R2, Lm2, Pxy
real rowvector V, Bt, V2, Bt2
transmorphic view_dv, view_iv
/*begin processing*/
st_view(view_dv, ., tokens(dvlist), st_varindex(touse))
st_view(view_iv, ., tokens(ivlist), st_varindex(touse))
Y = correlation(view_dv) //obtain DV correlations
X = correlation(view_iv) //obtain IV correlations
L = R = X //set-up for svd(); IV side
L2 = R2 = Y //set-up for svd(); DV side
V = J(1, cols(X), .) //placeholder for eigenvalues; IV side
V2 = J(1, cols(Y), .) //placeholder for eigenvalues; DV side
svd(X, L, V, R) //conduct singular value decomposition; IV side
svd(Y, L2, V2, R2) //conduct singular value decomposition; DV side
Lm = (L*diag(sqrt(V))*R) //process orthogonalized IVs
Lm2 = (L2*diag(sqrt(V2))*R2) //process orthogonalized DVs
Pxy = correlation((view_iv, view_dv)) //correlation between original IVs and DVs
Pxy = Pxy[rows(X)+1..rows(Pxy), 1..cols(X)] //take only IV-DV correlations
Bt2 = Pxy'*invsym(Lm2) //obtain adjusted DV interrelations
Bt = invsym(Lm)*Bt2 //obtain adjusted regression weights
Bt = Bt:^2 //square values of regression weights
Lm = Lm:^2 //square values of orthogonalized predictors
st_matrix("r(domwgts)", mean((Lm*Bt)')) //produce proportion of variance explained and put into Stata
st_numscalar("r(fs)", sum(mean((Lm*Bt)'))) //sum relative weights to obtain R2
}
end
/* programming notes and history
- mvdom version 1.0 - date - Jan 15, 2014
Basic version
-----
- mvdom version 1.1 - date - March 11, 2015
//notable changes\\
- added version statement (12.1)
- added the Pxy metric
- added the epsilon-based function
- changed canon to _canon; canon had odd behavior when called from mvdom
---
mvdom version 1.2.0 - August 14, 2023
- minimum versions 15 consistent with base -domin-
- 'if' statement optional to accommodate domin v. 3.5.0
- cleaned up returned values - consistent with reported
- added e(sample) needed for domin v. 3.5.0
- restructured eps_ri_mv()
- removed noconstant option - not consistently applied, not sure why one would ever use it
- unhide epsilon - note that it produces Pxy as a metric
- use st_view as opposed to st_data for efficiency
// 1.2.1 - December 20, 2024
- version to 16 consistent with base -domin-
** mi-able?