*! version 1.8.9 PR 21Dec98. program define mfracpol, eclass version 6 if "`1'" == "" | "`1'"=="," { if e(fp_cmd2)!="mfracpol" { error 301 } _fracrep "fractional polynomial" " df " "Powers" exit } gettoken cmd 0 : 0 frac_chk `cmd' if `s(bad)' { di in red "invalid or unrecognised command, `cmd'" exit 198 } *!! st commands disabled /* if $S_7 { _fracst m `cmd' `*' exit } */ /* dist=0 (normal), 1 (binomial), 2 (poisson), 3 (cox), 4 (glm), 5 (xtgee), 6(ereg/weibull). */ local dist `s(dist)' local glm `s(isglm)' local qreg `s(isqreg)' local xtgee `s(isxtgee)' local normal `s(isnorm)' syntax varlist(min=2) [if] [in] [aw fw pw iw] [, /* */ ADDpowers(str) ADJust(str) ALpha(str) ALL CATzero(str) DF(str) /* */ DFMax(int 4) CYCles(int 5) DEAD(str) FIXpowers(str) FP01 noSCAling /* */ POwers(str) SELect(str) SEQuential XOrder(str) XPowers(str) ZERo(str) * ] frac_cox "`dead'" `dist' /* Process options */ local regopt `options' if "`dead'"!="" { local regopt "`regopt' dead(`dead')" } if "`powers'"=="" { local powers "-2,-1,-.5,0,.5,1,2,3" } if "`addpowe'"!="" { local fpopt "`fpopt' addp(`addpowe')" } if "`fixpowe'"!="" { local fpopt "`fpopt' fixp(`fixpowe')" } /* Check for missing values in lhs, rhs and model vars. */ tempvar touse quietly { mark `touse' [`weight' `exp'] `if' `in' markout `touse' `varlist' `dead' frac_wgt "`exp'" `touse' "`weight'" local wgt `r(wgt)' /* [`weight'`exp'] */ count if `touse' local nobs = _result(1) } /* Rearrange order of variables in varlist */ if "`xorder'"=="" { local xorder "+" } _fracord `varlist' `wgt' if `touse', order(`xorder') `regopt' /* */ cmd(`cmd') `fp01' local varlist "$S_1" parse "`varlist'", parse(" ") local lhs `1' mac shift local rhs "`*'" /* Initialisation. Count RHS vars, store var names in n1, n2, ..., and store default df and powers. */ local nx 0 local l 0 while "`1'"!="" { local nx = `nx'+1 local alp`nx' .05 /* default FP selection level */ local h`nx' "`1'" /* names of H(xvars) */ local n`nx' "`1'" /* names of xvars */ local po`nx' 1 /* to be final power */ local sel`nx' 1 /* default var selection level */ /* Remove old I* variables for xvar */ frac_mun `n`nx'' purge /* Store original order of each RHS variable */ local i=`nx'+1 local o`nx' ${S_`i'} mac shift } /* Adjustment */ if "`adjust'"!="" { frac_dis "`adjust'" adjust "`rhs'" local j 0 while `j'<`nx' { local j=`j'+1 local adj`j' ${S_`j'} } } /* Set up degrees of freedom for each variable */ if "`df'"!="" { frac_dis "`df'" df 1 . "`rhs'" local j 0 while `j'<`nx' { local j=`j'+1 if "${S_`j'}"!="" { local df`j' ${S_`j'} } } } /* Assign default df for vars not so far accounted for. Award 1 df if 2-3 distinct values, 2 df for 4-5 values, dfmax df for >=6 values. Adjustment: adjust by default if >2 distinct values, not otherwise. */ local j 1 while `j'<=`nx' { qui inspect `n`j'' if `touse' local uniq=_result(7) if "`df`j''"=="" { local n `n`j'' if `uniq'==1 { noi di in red "`n`j'' does not vary." exit 198 } else if `uniq'<=3 { local df`j' 1 } else if `uniq'<=5 { local df`j' 2 } else local df`j' `dfmax' } if "`adj`j''"!="" { if "`adj`j''"=="no" { local adj`j' } else local adj`j' adjust(`adj') } else if `uniq'>2 { local adj`j' adjust(mean) } local j=`j'+1 } /* Set up FP selection level (alpha) for each variable */ if "`alpha'"!="" { frac_dis "`alpha'" alpha 0 1 "`rhs'" local j 0 while `j'<`nx' { local j=`j'+1 if "${S_`j'}"!="" { local alp`j' ${S_`j'} } } } /* Set up selection level for each variable */ if "`select'"!="" { frac_dis "`select'" select 0 1 "`rhs'" local j 0 while `j'<`nx' { local j=`j'+1 if "${S_`j'}"!="" { local sel`j' ${S_`j'} } } } /* Powers for individual FP terms. */ if "`xpowers'"!="" { parse "`xpowers'", parse(",") local ncl 0 /* # of comma-delimited clusters */ while "`1'"!="" { if "`1'"=="," { mac shift } local ncl=`ncl'+1 local clust`ncl' "`1'" mac shift } if "`clust`ncl''"=="" { local ncl=`ncl'-1 } if `ncl'>`nx' { di in red "too many sets of powers() values specified" exit 198 } /* Disentangle each varlist:string cluster */ local i 1 while `i'<=`ncl' { parse "`clust`i''", parse("=:") local powk `3' unabbrev "`1'" parse "$S_1", parse(" ") while "`1'"!="" { frac_in `1' "`rhs'" if `df`s(k)''<2 { di in red "xpowers invalid for `1', it has 1 df" exit 198 } local xpow$S_1 `powk' mac shift } local i = `i'+1 } } /* Vars with zero option */ if "`zero'"!="" { unabbrev "`zero'" parse "$S_1", parse(" ") while "`1'"!="" { frac_in `1' "`rhs'" local zero`s(k)' "zero" mac shift } } /* Vars with catzero option */ if "`catzero'"!="" { unabbrev "`catzero'" parse "$S_1", parse(" ") while "`1'"!="" { frac_in `1' "`rhs'" local catz`s(k)' "catzero" local zero`s(k)' "zero" /* catzero implies zero */ mac shift } } /* Reserve names for H(predictors) by creating a dummy variable for each predictor which potentially needs transformation. */ local i 1 while `i'<=`nx' { if `df`i''>1 | "`zero`i''`catz`i''"!="" { frac_mun `n`i'' local stub`i' `s(name)' qui gen byte `stub`i''_1=. } local i=`i'+1 } /* Build FP model. */ local it 0 local initial 1 local stable 0 /* convergence flag */ while !`stable' & `it'<=`cycles' { local it = `it'+1 local pwrs local rhs1 local stable 1 /* later changed to 0 if any power or status changes */ local lastch 0 /* becomes index of last var which changed status */ local i 1 while `i'<=`nx' { local ni "`n`i''" local dfi "df(`df`i'')" /* Build up RHS2 from the i+1th var to the end */ local rhs2 local j `i' while `j'<`nx' { local j = `j'+1 local rhs2 `rhs2' `h`j'' } if `initial' { if "`rhs2'"!="" { local fixed "base(`rhs2')" } else local fixed qui frac_sel `cmd' `lhs' `ni' `wgt' if `touse', /* */ df(1) `fixed' select(1) `regopt' di in gr _n /* */ "Deviance for model with all terms " /* */ "untransformed = " in ye %9.3f r(dev) in gr ", " /* */ in ye `nobs' in gr " observations" } if "`rhs1'`rhs2'"!="" { local fixed "base(`rhs1' `rhs2')" } else local fixed /* Vars with df(1) are straight-line */ local pvalopt "alpha(`alp`i'') select(`sel`i'')" if `i'==1 { di } if `df`i''==1 { local pw local fpo } else { if "`xpow`i''"!="" { local pw "powers(`xpow`i'')" } else local pw "powers(`powers')" local fpo `fpopt' } if `df`i''==1 & `sel`i''==1 { /* var is included anyway */ local rhs1 `rhs1' `h`i'' di in gr "[`ni' included with 1 df in model]" _n } else { if "`stub`i''"!="" { local n name(`stub`i'') } else local n frac_sel `cmd' `lhs' `ni' `wgt' if `touse', `dfi' /* */ `pw' `zero`i'' `catz`i'' `fixed' `h' /* */ `fpo' `regopt' `pvalopt' `n' `sequent' local h`i' `r(rhs)' local dev=r(dev) local p `r(pwrs)' if "`p'"!="`po`i''" { if `nx'>1 { local stable 0 } local po`i' "`p'" local lastch `i' } if "`pwrs'"=="" { local pwrs "`p'" } else local pwrs "`pwrs',`p'" if "`h`i''"!="" { local rhs1 `rhs1' `h`i'' } } if `initial' { local h "nohead" local initial 0 } local i = `i'+1 } if `lastch'==1 { local stable 1 } /* 1 change only, at i=1 */ if !`stable' { di in gr _dup(70) "-" _n "Cycle " in ye `it' in gr /* */ ": deviance = " in ye %9.3f `dev' _n in gr _dup(70) "-" } } if `nx'>1 { local s if `it'!=1 { local s "s" } if !`stable' { di _n in gr "No convergence" _cont } else di _n in gr /* */ "Fractional polynomial fitting algorithm converged" _cont di in gr " after " in ye `it' in gr " cycle`s'." } if `stable' { di _n in gr "Transformations of covariates:" _n } /* Remove variables left behind by _frac154 */ local i 1 while `i'<=`nx' { if "`stub`i''"!="" { cap drop `stub`i''* } local i=`i'+1 } /* Store results */ local finalvl /* predictors in final model */ local i 1 while `i'<=`nx' { local o `o`i'' local p `po`o'' local x `n`o'' local z `zero`o'' local c `catz`o'' local a `adj`o'' /* Create FP vars as necessary, with new unique names */ if trim("`p'")!="." { if trim("`p'")=="1" & "`z'`c'`a'"=="" { local namex `x' } else { frac_mun `x' local vn `s(name)' fracgen `x' `p' if e(sample), `all' `scaling' /* */ `z' `c' `a' name(`vn') local namex `r(names)' } local finalvl `finalvl' `namex' } /* Save stuff for _fracrep to use for display purposes */ local npows: word count `p' local fdf=2*`npows'-(`npows'==1&"`p'"=="1")+("`c'"!="") local idf=`df`o''+("`c'"!="") global MFPfd`i' `fdf' /* final degrees of freedom */ global MFPdf`i' `idf' /* initial degrees of freedom */ global MFPal`i' `alp`o'' /* FP selection level */ global MFPse`i' `sel`o'' /* var selection levl */ local i=`i'+1 } /* Estimate final (conditionally linear) model. */ qui `cmd' `lhs' `finalvl' `wgt' if e(sample), `regopt' /* Store results */ global S_1 `finalvl' global S_2 `dev' local i 0 while `i'<`nx' { local i=`i'+1 local o `o`i'' global S_E_x`i' `n`o'' /* name of ith predictor in user order */ est local fp_x`i' `n`o'' global S_E_k`i' `po`o'' /* "powers" for ith predictor */ est local fp_k`i' `po`o'' if "`catz`o''"!="" { global S_E_c`i' 1 est local fp_c`i' 1 } } global S_E_dist `dist' est scalar fp_dist=`dist' global S_E_wgt `weight' est local fp_wgt `weight' global S_E_exp "`exp'" est local fp_exp `exp' global S_E_depv `lhs' est local fp_depv `lhs' global S_E_dev `dev' est scalar fp_dev=`dev' global S_E_rhs /* deliberately blank for consistency with fracpoly */ est local fp_rhs global S_E_opts `regopt' est local fp_opts `regopt' global S_E_fvl `finalvl' est local fp_fvl `finalvl' global S_E_nx `nx' est scalar fp_nx=`nx' global S_E_t1t "Fractional Polynomial" est local fp_t1t "Fractional Polynomial" _fracrep "fractional polynomial" " df " "Powers" global S_E_fp "fracpoly" est local fp_cmd "fracpoly" global S_E_fp2 "mfracpol" est local fp_cmd2 "mfracpol" end