Skip to content

Instantly share code, notes, and snippets.

@doctorkid
Created June 12, 2018 18:14
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save doctorkid/b6db751ff947b825265f70638ce18b61 to your computer and use it in GitHub Desktop.
Save doctorkid/b6db751ff947b825265f70638ce18b61 to your computer and use it in GitHub Desktop.
*BC 2nd
clear
input time n_male n_female p_bw_male p_bw_female p_bd_male p_bd_female p_bn_male p_bn_female p_all
9 34 38 0.08 0.12 0 0 0 0 0
12 46 50 0.32 0.39 0 0 0 0 0
18 85 79 0.57 0.65 0.10 0.16 0.02 0.06 0.028
24 152 150 0.73 0.83 0.17 0.22 0.07 0.09 0.054
36 162 161 0.96 0.98 0.18 0.18 0.16 0.22 0.109
48 163 161 0.98 0.97 0.89 0.89 0.86 0.86 0.766
60 143 146 0.98 0.99 0.92 0.92 0.86 0.89 0.855
72 152 156 0.99 0.99 0.97 0.97 0.90 0.94 0.913
end
gen pct_bw_male = 100*p_bw_male
gen pct_bw_female = 100*p_bw_female
gen pct_bd_male = 100*p_bd_male
gen pct_bd_female = 100*p_bd_female +1
gen pct_bn_male = 100*p_bn_male
gen pct_bn_female = 100*p_bn_female
gen pct_all = 100*p_all
*Total
gen N = n_male + n_female
gen N_bw = n_male*p_bw_male + n_female*p_bw_female
gen Pct_bw = N_bw/N*100
gen N_bd = n_male*p_bd_male + n_female*p_bd_female
gen Pct_bd = N_bd/N*100
gen N_bn = n_male*p_bn_male + n_female*p_bn_female
gen Pct_bn = N_bn/N*100
graph twoway (line Pct_bw time, lcolor(gray) lwidth(thin) ) (line Pct_bd time, color(orange) lwidth(thin) ) (line Pct_bn time, color(purple) lwidth(thin) ) (line pct_all time, color(black) lwidth(thin)) ///
(scatter Pct_bw time, color(gray) msize(small)) (scatter Pct_bd time, color(orange) msize(small)) (scatter Pct_bn time, color(purple) msize(small)) (scatter pct_all time, color(black) msize(small)) ///
, legend(off) ///
xtitle(月齢) ytitle(排便・排尿コントロール完了 (%)) ///
xlabel(6(6)72) xline(12 24 36 48 60 72, lcolor(gray) lpattern(dot)) ///
graphregion( color(white) ) plotregion( fcolor(white) )
gen p_bw = Pct_bw/100
gen p_bd = Pct_bd/100
gen p_bn = Pct_bn/100
gen ln_odds_bw = ln(p_bw)/(1-ln(p_bw))
gen ln_odds_bd = ln(p_bd)/(1-ln(p_bd))
gen ln_odds_bn = ln(p_bn)/(1-ln(p_bn))
gen ln_odds_all = ln(p_all)/(1-ln(p_all))
gen ln_time = ln(time)
glm p_bw c.time ln_time [fweight=N], link(logit) family(binomial)
glm p_bd c.time ln_time [fweight=N], link(logit) family(binomial)
glm p_bn c.time ln_time [fweight=N], link(logit) family(binomial)
glm p_all c.time ln_time [fweight=N], link(logit) family(binomial)
clear
set obs 63
gen time = int((_n-1)) +9
gen ln_time = ln(time)
gen ln_odds_bw = -9.654204-.0097237*time + 3.570269*ln_time
gen ln_odds_bd = -2.607611+.1533049 *time -.9545885*ln_time
gen ln_odds_bn = -19.97697-.0021136 *time+ 5.411243*ln_time
gen ln_odds_all = -21.90521-.0017088 *time+ 5.816246*ln_time
gen p_bw = exp(ln_odds_bw)/(1+exp(ln_odds_bw))*100
gen p_bd = exp(ln_odds_bd)/(1+exp(ln_odds_bd))*100
gen p_bn = exp(ln_odds_bn)/(1+exp(ln_odds_bn))*100
gen p_all = exp(ln_odds_all)/(1+exp(ln_odds_all))*100
graph twoway (line p_bw time, lcolor(gray) lwidth(thin) ) (line p_bd time, color(orange) lwidth(thin) ) (line p_bn time, color(purple) lwidth(thin) ) (line p_all time, color(black) lwidth(thin)) ///
(scatter p_bw time, color(gray) msize(vsmall)) (scatter p_bd time, color(orange) msize(vsmall)) (scatter p_bn time, color(purple) msize(vsmall)) (scatter p_all time, color(black) msize(vsmall)) ///
, legend(off) ///
xtitle(月齢) ytitle(排便・排尿コントロール完了 (%)) ///
xlabel(6(6)72) xline(12 24 36 48 60 72, lcolor(gray) lpattern(dot)) ///
graphregion( color(white) ) plotregion( fcolor(white) )
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment