* version 4.3, DEC2013 * Artur Pokropek & Maciej Jakubowski * regression with PIAAC data cap program drop piaacreg program define piaacreg, rclass syntax [varlist(default=none)] [if] [in] /// , save(string) [ weight(varlist max=1) rep(varlist) /// countryid(varname) vemethodn(varname numeric) pvdep(string) pvindep1(string) pvindep2(string) pvindep3(string) /// over(string) fast round(integer 2) cmd(string) cmdops(string) r2(string) cons or] version 10.0 tempvar probka cntid jkvar proflevel1 tempname b se over_values local oldover="`over'" local longlist=0 local overpv=0 if "`over'"=="pvlit" | "`over'"=="pvnum" | "`over'"=="PVLIT" | "`over'"=="PVNUM" { qui recode `over'1 (.=.) (0/175.9999=0) (176/225.9999=1) (226/275.9999=2) (276/325.9999=3) (326/375.9999=4) (376/999=5), gen(`proflevel1') local overpv=1 local over="`proflevel1'" } else if "`over'"=="pvpsl" | "`over'"=="PVPLS" { qui recode pvpsl1 (.=.) (0/240.9999=0) (241/290.9999=1) (291/340.9999=2) (341/999=3), gen(`proflevel1') local overpv=1 local over="`proflevel1'" } if "`cmd'"=="" local cmd="reg" if "`r2'"=="" { if "`cmd'"=="reg" local r2="r2" else if "`cmd'"=="logit" | "`cmd'"=="logistic" | "`cmd'"=="poisson" | "`cmd'"=="ologit" local r2="r2_p" } if "`pvdep'"=="" { tokenize `varlist' local outcome `1' macro shift local varlist `*' local pv=0 } else { local outcome `pvdep'1 local pv=1 } if "`cmd'"=="logit" | "`cmd'"=="logistic" { sum `outcome', meanonly if r(min)!=0 | r(mean)<=0 { di as error "outcome has to vary and contain 0 and positive values only" } } if "`cmd'"=="poisson" { sum `outcome', meanonly if r(mean)<=0 { di as error "outcome has to vary and contain nonnegative values only" } } if "`pvindep1'"!="" | "`pvindep2'"!="" | "`pvindep3'"!="" | `overpv'!=0 local pv=1 if "`countryid'"!="" { capture confirm string variable `countryid' if !_rc { qui gen `cntid'=`countryid' qui levelsof `cntid' `if' `in', local(cnt) } else { qui tostring(`countryid'), gen(`cntid') qui levelsof `cntid' `if' `in', local(cnt) } } else { cap confirm variable cntryid if !_rc { qui tostring(cntryid), gen(`cntid') qui local cnt="36 40 124 203 208 233 246 250 276 372 380 392 410 528 578 616 703 724 752 840 SUB 56 826 OECD 196 643" local longlist=1 } else { cap confirm variable CNTRYID if !_rc { qui tostring(CNTRYID), gen(`cntid') qui local cnt="36 40 124 203 208 233 246 250 276 372 380 392 410 528 578 616 703 724 752 840 SUB 56 826 OECD 196 643" local longlist=1 } else { di as error "Your dataset does not contain variable cntryid or CNTRYID. Please specify option countryid()" exit 100 } } } if "`weight'"=="" { cap confirm variable spfwt0 if !_rc local weight="spfwt0" else { cap confirm variable SPFWT0 if !_rc local weight="SPFWT0" else { di as error "Your dataset does not contain main survey weight spfwt0 or SPFWT0. Please specify weight() option" exit 100 } } } local ok_small=0 local ok_big=0 if "`rep'"=="" { forvalues i=1(1)80 { cap confirm variable spfwt`i' if !_rc local ok_small=`ok_small'+1 cap confirm variable SPFWT`i' if !_rc local ok_big=`ok_big'+1 } if `ok_small'==80 local rep="spfwt1-spfwt80" else if `ok_big'==80 local rep="SPFWT1-SPFWT80" else { di as error "Your dataset does not contain all replicate weights spfwt1-spfwt80 or SPFWT1-SPFWT80. Please specify rep() option" exit 100 } } cap confirm variable vemethodn if !_rc { qui gen `jkvar'=vemethodn } else { cap confirm variable VEMETHODN if !_rc { qui gen `jkvar'=VEMETHODN } else { if "`vemethodn'"=="" { di as error "Your dataset does not contain variable vemethodn or VEMETHODN. Please specify option vemethodn()" } else gen `jkvar'=`vemethodn' } } if "`if'"=="" local if=" if " else local if=`"`if' "'+" & " tempvar test qui gen `test'=1 foreach var of local varlist { capture confirm numeric variable `var' if !_rc qui replace `test'=`test'*`var' else di as error "`var' is not a numeric variable" } if "`pvdep'"=="" { capture confirm numeric variable `outcome' if !_rc qui replace `test'=`test'*`outcome' else di as error "`outcome' is not a numeric variable" } else { qui replace `test'=`test'*`pvdep'1 local outcome="" } if "`pvdep'"!="" qui replace `test'=`test'*`pvdep'1 if "`pvindep1'"!="" qui replace `test'=`test'*`pvindep1'1 if "`pvindep2'"!="" qui replace `test'=`test'*`pvindep2'1 if "`pvindep3'"!="" qui replace `test'=`test'*`pvindep3'1 if `overpv'==1 qui replace `test'=`test'*`oldover'1 local n_vars: word count `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' local n_cnt: word count `cnt' if `longlist'==1 local n_cnt=`n_cnt'-2 if "`over'"!="" { tab `over', nofreq matrow(over_values) local n_over=r(r) } else { local over="1" matrix over_values=(1,1) local n_over=1 } tempname tempfile file open `tempfile' using "`save'.html", write replace file write `tempfile' `""' " " * TABLE HEADER - OVER CATEGORIES - if over() specified if `n_over'>1 { local col_span=2*`n_vars'+1 file write `tempfile' `""' forvalues i=1(1)`n_over' { local numerek=over_values[`i',1] if `overpv'!=1 { local labelka: label (`over') `numerek' file write `tempfile' `""' } else { if `numerek'==0 { file write `tempfile' `""' } else { file write `tempfile' `""' } } } } * SECOND ROW - predictors file write `tempfile' `"
Over categories:`labelka'Below Level 1Level `numerek'
"' forvalues i=1(1)`n_over' { foreach var in `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' { file write `tempfile' `""' } file write `tempfile' `""' } * third row file write `tempfile' "" forvalues i=1(1)`n_over' { foreach var in `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' { file write `tempfile' `"" `"" local OECD_`var'_`i'=0 local OECD_SE_`var'_`i'=0 } file write `tempfile' `"" local OECD_r2_`i'=0 } local decimal=0.1^`round' forvalues i=1(1)`n_over' { local noc_`i'=0 * creating output matrices mat b_over`i'=J(`n_cnt',`n_vars',.) mat se_over`i'=J(`n_cnt',`n_vars',.) mat r2_over`i'=J(`n_cnt',1,.) } local OECD_saved=0 local pvlist="`outcome'" if "`pvindep1'"!="" local pvlist="`pvlist' `pvindep1'*" if "`pvindep2'"!="" local pvlist="`pvlist' `pvindep2'*" if "`pvindep3'"!="" local pvlist="`pvlist' `pvindep3'*" if "`pvdep'"!="" local pvlist="`pvlist' `pvdep'*" if `overpv'==1 local pvlist="`pvlist' `oldover'*" local country_list="" local cnt_i=0 foreach l of local cnt { if "`l'"=="OECD" | "`l'"=="SUB" { if "`l'"=="SUB" file write `tempfile' `""' else { local OECD_saved=1 file write `tempfile' `""' forvalues i=1(1)`n_over' { foreach var in `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' { local mean=string(round(`OECD_`var'_`i''/`noc_`i'',`decimal'),"%20.`round'f") local se=string(round(sqrt(`OECD_SE_`var'_`i'')/`noc_`i'',`decimal'),"%20.`round'f") file write `tempfile' `" "' } local mean=string(round(`OECD_r2_`i''/`noc_`i'',`decimal'),"%20.`round'f") file write `tempfile' `" "' } file write `tempfile' `""' } } else { _cnt `l' local name=r(name) local country_list="`country_list'`l' " local cnt_i=`cnt_i'+1 di "" di as result "`name'" _continue file write `tempfile' "" sum `jkvar' if `cntid'=="`l'" , meanonly if `r(N)'>0 { local jk=1 sum `jkvar' if `cntid'=="`l'" , meanonly cap local jk=`r(mean)' if `jk'<1 | `jk'>2 { di as error "Please specify values 1 or 2 in `jkvar'" exit 100 } forvalues i=1(1)`n_over' { if 1!=`n_over' di as text " `i'" _c sum `test' if `cntid'=="`l'" & `over'==over_values[`i',1], meanonly if r(N)>2 { return clear cap drop `probka' if `overpv'==0 qui gen `probka'=1 `if' `cntid'=="`l'" & `over'==over_values[`i',1] `in' else qui gen `probka'=1 `if' `cntid'=="`l'" `in' local temp=over_values[`i',1] preserve qui keep if `probka'==1 keep `pvlist' `varlist' `probka' `weight' `rep' if `pv'==1 | `overpv'==1 { if "`fast'"=="" _piaacregpv `outcome' `varlist', overval(`temp') overpv(`oldover') if(`probka') pvdep(`pvdep') pvindep1(`pvindep1') pvindep2(`pvindep2') pvindep3(`pvindep3') weight(`weight') rep(`rep') jk(`jk') cmd(`cmd') cmdops(`cmdops') r2(`r2') `cons' `or' else _piaacregfastpv `outcome' `varlist', overval(`temp') overpv(`oldover') if(`probka') pvdep(`pvdep') pvindep1(`pvindep1') pvindep2(`pvindep2') pvindep3(`pvindep3') weight(`weight') cmd(`cmd') cmdops(`cmdops') r2(`r2') `cons' `or' } else { if "`fast'"=="" _piaacreg `outcome' `varlist', if(`probka') weight(`weight') rep(`rep') jk(`jk') cmd(`cmd') cmdops(`cmdops') r2(`r2') `cons' `or' else _piaacregfast `outcome' `varlist', if(`probka') weight(`weight') cmd(`cmd') cmdops(`cmdops') r2(`r2') `cons' `or' } matrix b=r(coef) matrix se=r(se) local noc_`i'=`noc_`i''+1 mat temp=b' mat temp= temp[1,1..`n_vars'] mat b_over`i'[`cnt_i',1]=temp mat temp=se' mat temp= temp[1,1..`n_vars'] mat se_over`i'[`cnt_i',1]=temp local row=1 foreach var in `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' { local b=b[`row',1] local se=se[`row',1] local row=`row'+1 if `b'==0 & `se'==0 file write `tempfile' `""' `""' else { local mean=string(round(`b',`decimal'),"%20.`round'f") local se=string(round(`se',`decimal'),"%20.`round'f") file write `tempfile' `""' `""' local OECD_`var'_`i'=`OECD_`var'_`i''+`mean' local OECD_SE_`var'_`i'=`OECD_SE_`var'_`i''+`se'^2 } } local b=b[`row',1] mat r2_over`i'[`cnt_i',1]=`b' local mean=string(round(`b',`decimal'),"%20.`round'f") file write `tempfile' `""' local OECD_r2_`i'=`OECD_r2_`i''+`mean' restore } else { foreach var in `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' { file write `tempfile' `""' `""' } file write `tempfile' `""' } } } } } if `OECD_saved'==0 { file write `tempfile' `""' forvalues i=1(1)`n_over' { foreach var in `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' { local mean=string(round(`OECD_`var'_`i''/`noc_`i'',`decimal'),"%20.`round'f") local se=string(round(sqrt(`OECD_SE_`var'_`i'')/`noc_`i'',`decimal'),"%20.`round'f") file write `tempfile' `" "' } local mean=string(round(`OECD_r2_`i''/`noc_`i'',`decimal'),"%20.`round'f") file write `tempfile' `" "' } } file write `tempfile' _n "
Variable`var'`r2'
" "Country" ""' "Coef." ""' "S.E." ""' "" "
Sub-national entities
OECD Average `mean' `se' `mean'
Partners
" "`name'" " 0 0 `mean' `se' `mean' . . .
Average `mean' `se' `mean'
" file close `tempfile' di "" di "Results saved in the `save'.html file" forvalues i=1(1)`n_over' { matrix colnames b_over`i' = `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' matrix colnames se_over`i' = `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' matrix rownames b_over`i' = `country_list' matrix rownames se_over`i' = `country_list' matrix rownames r2_over`i' = `country_list' matrix colnames r2_over`i' = `r2' } if `n_over'==1 { return matrix b=b_over1 return matrix se=se_over1 return matrix r2=r2_over1 } else { forvalues i=1(1)`n_over' { return matrix b_over`i'=b_over`i' return matrix se_over`i'=se_over`i' return matrix r2_over`i'=r2_over`i' } } end *********** cap program drop _cnt program define _cnt, rclass args l local name="`l'" if "`l'"=="36" local name="Australia" if "`l'"=="40" local name="Austria" if "`l'"=="124" local name="Canada" if "`l'"=="203" local name="Czech Republic" if "`l'"=="276" local name="Germany" if "`l'"=="208" local name="Denmark" if "`l'"=="724" local name="Spain" if "`l'"=="233" local name="Estonia" if "`l'"=="246" local name="Finland" if "`l'"=="250" local name="France" if "`l'"=="826" local name="England/N. Ireland (UK)" if "`l'"=="372" local name="Ireland" if "`l'"=="380" local name="Italy" if "`l'"=="392" local name="Japan" if "`l'"=="410" local name="Korea" if "`l'"=="528" local name="Netherlands" if "`l'"=="578" local name="Norway" if "`l'"=="616" local name="Poland" if "`l'"=="643" local name="Russian Federation" if "`l'"=="703" local name="Slovak Republic" if "`l'"=="752" local name="Sweden" if "`l'"=="840" local name="United States" if "`l'"=="56" local name="Flanders (Belgium)" if "`l'"=="196" local name="Cyprus" return local name "`name'" end **************************************************** cap program drop _piaacreg program define _piaacreg, rclass syntax [varlist], weight(varlist max=1) rep(varlist) jk(integer) /// cmd(string) if(varname numeric) [cmdops(string) r2(string) cons or] tokenize `varlist' local outcome `1' macro shift local varlist `*' * SETUP local num_rep : word count `rep' local n_vars: word count `varlist' if "`cons'"!="" { local n_vars=`n_vars'+1 local stala="_cons" } matrix coef = J(`n_vars'+1,1,.) matrix rownames coef =`varlist' `cons' `r2' matrix se= J(`n_vars',1,.) matrix rownames se =`varlist' `cons' qui `cmd' `outcome' `varlist' [pw=`weight'] if `if'==1, `cmdops' foreach var in `varlist' `stala' { local fulstat_`var'=_b[`var'] local sumofsq_`var'=0 } local fulstat_r2=e(`r2') foreach replication of local rep { qui `cmd' `outcome' `varlist' [pw=`replication'] if `if'==1, `cmdops' foreach var in `varlist' `stala' { local sq_`var' = (_b[`var'] - `fulstat_`var'')^2 local sumofsq_`var'=`sumofsq_`var''+ `sq_`var'' } } foreach var in `varlist' `stala' { if `jk'==1 local sdr_`var'=sqrt(((`num_rep'-1)/`num_rep')*`sumofsq_`var'') else local sdr_`var'=sqrt(1*`sumofsq_`var'') } local row=1 foreach var in `varlist' `stala' { if "`or'"=="" { matrix coef[`row',1] = `fulstat_`var'' matrix se[`row',1] = `sdr_`var'' } else { matrix coef[`row',1] = exp(`fulstat_`var'') matrix se[`row',1] = exp(`fulstat_`var'')*`sdr_`var'' } local row=`row'+1 } matrix coef[`row',1] = `fulstat_r2' matrix colnames coef = b matrix colnames se = se return matrix coef =coef return matrix se =se end ************************************* cap program drop _piaacregfast program define _piaacregfast, rclass syntax [varlist], weight(varlist max=1) /// cmd(string) if(varname numeric) [cmdops(string) r2(string) cons or] tokenize `varlist' local outcome `1' macro shift local varlist `*' * SETUP local n_vars: word count `varlist' if "`cons'"!="" { local n_vars=`n_vars'+1 local stala="_cons" } matrix coef = J(`n_vars'+1,1,.) matrix rownames coef =`varlist' `cons' `r2' matrix se= J(`n_vars',1,.) matrix rownames se =`varlist' `cons' qui `cmd' `outcome' `varlist' [pw=`weight'] if `if'==1, `cmdops' local row=1 foreach var in `varlist' `stala' { if "`or'"=="" { matrix coef[`row',1] = _b[`var'] matrix se[`row',1] = _se[`var'] } else { matrix coef[`row',1] = exp(_b[`var']) matrix se[`row',1] = exp(_b[`var'])*_se[`var'] } local row=`row'+1 } matrix coef[`row',1] = e(`r2') matrix colnames coef = b matrix colnames se = se return matrix coef =coef return matrix se =se end ************** cap program drop _piaacregpv program define _piaacregpv, rclass syntax [varlist(numeric default=none)], [overpv(string) overval(string) pvdep(string) pvindep1(string) pvindep2(string) pvindep3(string)] /// weight(varlist max=1) rep(varlist) jk(integer) cmd(string) if(varname numeric) [cmdops(string) r2(string) cons or] if "`overpv'"=="pvlit" | "`overpv'"=="pvnum" | "`overpv'"=="PVLIT" | "`overpv'"=="PVNUM" { tempvar proflevel1 proflevel2 proflevel3 proflevel4 proflevel5 proflevel6 proflevel7 proflevel8 proflevel9 proflevel10 forvalues i=1(1)10 { qui recode `overpv'`i' (.=.) (0/175.9999=0) (176/225.9999=1) (226/275.9999=2) (276/325.9999=3) (326/375.9999=4) (376/999=5), gen(`proflevel`i'') } local overpv=1 } else if "`overpv'"=="pvpsl" | "`overpv'"=="PVPSL" { tempvar proflevel1 proflevel2 proflevel3 proflevel4 proflevel5 proflevel6 proflevel7 proflevel8 proflevel9 proflevel10 forvalues i=1(1)10 { qui recode pvpsl`i' (.=.) (0/240.9999=0) (241/290.9999=1) (291/340.9999=2) (341/999=3), gen(`proflevel`i'') } local overpv=1 } else local overpv=0 tempvar probkaif if "`pvdep'"=="" { tokenize `varlist' local outcome `1' macro shift local varlist `*' local pv=0 } else local pv=1 local num_rep : word count `rep' local n_vars: word count `varlist' `pvindep1' `pvindep2' `pvindep3' if "`cons'"!="" { local n_vars=`n_vars'+1 local stala="_cons" } matrix coef = J(`n_vars'+1,10,.) matrix rownames coef =`varlist' `pvindep1' `pvindep2' `pvindep3' `cons' `r2' matrix se= J(`n_vars',10,.) matrix rownames se =`varlist' `pvindep1' `pvindep2' `pvindep3' `cons' forvalues ipv=1(1)10 { di as text "." _c if "`pvindep1'"!="" local pvindep="`pvindep1'`ipv'" if "`pvindep2'"!="" local pvindep="`pvindep' `pvindep2'`ipv'" if "`pvindep3'"!="" local pvindep="`pvindep' `pvindep3'`ipv'" if `pv'!=0 local outcome="`pvdep'`ipv'" cap drop `probkaif' if `overpv'==1 qui gen `probkaif'=1 if `proflevel`ipv''==`overval' else local probkaif=1 qui `cmd' `outcome' `pvindep' `varlist' [pw=`weight'] if `probkaif'==1, `cmdops' foreach var in `pvindep' `varlist' `stala' { local fulstat_`var'=_b[`var'] local sumofsq_`var'=0 } local fulstat_r2=e(`r2') if `fulstat_r2'==. local fulstat_r2=0 foreach replication of local rep { qui `cmd' `outcome' `pvindep' `varlist' [pw=`replication'] if `probkaif'==1, `cmdops' foreach var in `pvindep' `varlist' `stala' { local sq_`var' = (_b[`var'] - `fulstat_`var'')^2 local sumofsq_`var'=`sumofsq_`var''+ `sq_`var'' } } foreach var in `pvindep' `varlist' `stala' { if `jk'==1 local rVAR_`var'=(((`num_rep'-1)/`num_rep')*`sumofsq_`var'') else local rVAR_`var'=(1*`sumofsq_`var'') } local row=1 foreach var in `pvindep' `varlist' `stala' { matrix coef[`row',`ipv'] = `fulstat_`var'' matrix se[`row',`ipv'] = `rVAR_`var'' local row=`row'+1 } matrix coef[`row',`ipv'] = `fulstat_r2' } mat coef=coef' mat se=se' mata: pvmat() if "`or'"!="" { local row=1 foreach var in `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' { matrix newcoef[`row',1] = exp(newcoef[`row',1]) matrix newse[`row',1] = newcoef[`row',1]*newse[`row',1] local row=`row'+1 } } matrix rownames newcoef = `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' `r2' matrix colnames newcoef = b matrix rownames newse = `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' matrix colnames newse = se return matrix coef =newcoef return matrix se =newse end ************** cap program drop _piaacregfastpv program define _piaacregfastpv, rclass syntax [varlist(numeric default=none)], [overpv(string) overval(string) pvdep(string) pvindep1(string) pvindep2(string) pvindep3(string)] /// weight(varlist max=1) cmd(string) if(varname numeric) [cmdops(string) r2(string) cons or] if "`overpv'"=="pvlit" | "`overpv'"=="pvnum" | "`overpv'"=="PVLIT" | "`overpv'"=="PVNUM" { tempvar proflevel1 proflevel2 proflevel3 proflevel4 proflevel5 proflevel6 proflevel7 proflevel8 proflevel9 proflevel10 forvalues i=1(1)10 { qui recode `overpv'`i' (.=.) (0/175.9999=0) (176/225.9999=1) (226/275.9999=2) (276/325.9999=3) (326/375.9999=4) (376/999=5), gen(`proflevel`i'') } local overpv=1 } else if "`overpv'"=="pvpsl" | "`overpv'"=="PVPSL" { tempvar proflevel1 proflevel2 proflevel3 proflevel4 proflevel5 proflevel6 proflevel7 proflevel8 proflevel9 proflevel10 forvalues i=1(1)10 { qui recode pvpsl`i' (.=.) (0/240.9999=0) (241/290.9999=1) (291/340.9999=2) (341/999=3), gen(`proflevel`i'') } local overpv=1 } else local overpv=0 tempvar probkaif if "`pvdep'"=="" { tokenize `varlist' local outcome `1' macro shift local varlist `*' local pv=0 } else local pv=1 local n_vars: word count `varlist' `pvindep1' `pvindep2' `pvindep3' if "`cons'"!="" { local n_vars=`n_vars'+1 local stala="_cons" } matrix coef = J(`n_vars'+1,10,.) matrix rownames coef =`varlist' `pvindep1' `pvindep2' `pvindep3' `cons' `r2' matrix se= J(`n_vars',10,.) matrix rownames se =`varlist' `pvindep1' `pvindep2' `pvindep3' `cons' forvalues ipv=1(1)10 { di as text "." _c if "`pvindep1'"!="" local pvindep="`pvindep1'`ipv'" if "`pvindep2'"!="" local pvindep="`pvindep' `pvindep2'`ipv'" if "`pvindep3'"!="" local pvindep="`pvindep' `pvindep3'`ipv'" if `pv'!=0 local outcome="`pvdep'`ipv'" cap drop `probkaif' if `overpv'==1 qui gen `probkaif'=1 if `proflevel`ipv''==`overval' else local probkaif=1 qui `cmd' `outcome' `pvindep' `varlist' [pw=`weight'] if `probkaif'==1, `cmdops' local row=1 foreach var in `pvindep' `varlist' `stala' { if "`or'"=="" { matrix coef[`row',`ipv'] = _b[`var'] matrix se[`row',`ipv'] = _se[`var'] } else { matrix coef[`row',`ipv'] = exp(_b[`var']) matrix se[`row',`ipv'] = exp(_b[`var'])*_se[`var'] } local row=`row'+1 } matrix coef[`row',`ipv'] = e(`r2') } mat coef=coef' mat se=se' mata: pvmat() matrix rownames newcoef = `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' `r2' matrix colnames newcoef = b matrix rownames newse = `pvindep1' `pvindep2' `pvindep3' `varlist' `cons' matrix colnames newse = se return matrix coef =newcoef return matrix se =newse end *************** version 10 clear mata mata: void pvmat() { pvs_org=st_matrix("coef") pvs=pvs_org[|1,1\10,cols(pvs_org)-1|] repVAR=mean(st_matrix("se")) pvdif=pvs:-mean(pvs) sumofsqpv=pvdif'pvdif VAR=1.1:*sumofsqpv:/9 pvVAR= diagonal(VAR)' totVAR=repVAR+pvVAR SDR=totVAR:^0.5 st_matrix("newse", SDR') st_matrix("newcoef", mean(pvs_org)') } end