. do "C:\Users\DEFAUL~1.SID\AppData\Local\Temp\STD00000000.tmp"

. * do-file for lecture 3b of VHM 802/812, Winter 2014
. version 13

. set more off

. cd "h:\vhm\vhm802\data_stata"
h:\vhm\vhm802\data_stata

.
. use mice.dta, clear

. regress dead dose /* to illustrate extrapolation below zero in regression */

Source | SS df MS Number of obs = 120

+
F( 1, 118) = 27.41
Model | 4.9616367 1 4.9616367 Prob > F = 0.0000
Residual | 21.3633633 118 .181045452 R-squared = 0.1885

+
Adj R-squared = 0.1816
Total | 26.325 119 .221218487 Root MSE = .42549

+

--------

dead | Coef. Std. Err. t P>|t| [95% Con f. Interval]

+

--------

dose | 2.91441 .5567137 5.24 0.000 1.811965 4.016855
_cons | -.2107171 .1094569 -1.93 0.057 -.4274715 .0060373

+

--------

.
. * logistic regression using logit command
. logit dead dose

Iteration 0: log likelihood = -75.669723
Iteration 1: log likelihood = -64.201632
Iteration 2: log likelihood = -63.944944
Iteration 3: log likelihood = -63.944713
Iteration 4: log likelihood = -63.944713

Logistic regression Number of obs = 120
LR chi2(1) = 23.45
Prob > chi2 = 0.0000
Log likelihood = -63.944713 Pseudo R2 = 0.1549

+

--------

dead | Coef. Std. Err. z P>|z| [95% Con f. Interval]

+

--------

dose | 14.6369 3.332533 4.39 0.000 8.10526 21.16855
_cons | -3.569738 .7053432 -5.06 0.000 -4.952185 -2.187291

+

--------

. predict linpred, xb /* predictions on logit scale */

. twoway (scatter linpred dose) (line linpred dose)

. lincom _cons+dose*0.1

( 1) .1*[dead]dose + [dead]_cons = 0

+

--------

dead | Coef. Std. Err. z P>|z| [95% Con f. Interval]

+

--------

(1) | -2.106048 .4015256 -5.25 0.000 -2.893023 -1.319072

+

--------

. predict phat, p /* predicted probabilities */

. twoway (scatter phat dose) (line phat dose)

. estat gof /* Pearson goodness-of-fit test */

Logistic model for dead, goodness-of-fit test

number of observations = 120
number of covariate patterns = 12
Pearson chi2(10) = 8.74
Prob > chi2 = 0.5567

.
. * analysis with dichotomous dose
. prtesti 60 32 60 7, count /* 2 independent proportions */

Two-sample test of proportions x: Number of obs = 60
y: Number of obs = 60

+

Variable | Mean Std. Err. z P>|z| [95% Conf. Interval]

+

x | .5333333 .0644061 .4070997 .659567
y | .1166667 .0414438 .0354382 .1978951

+

diff | .4166667 .0765881 .2665567 .5667766
| under Ho: .0855132 4.87 0.000

+

diff = prop(x) - prop(y) z = 4.8725
Ho: diff = 0

Ha: diff < 0 Ha: diff != 0 Ha: diff > 0
Pr(Z < z) = 1.0000 Pr(|Z| < |z|) = 0.0000 Pr(Z > z) = 0.0000

. csi 32 7 28 53, or /* relative risk and odds-ratio */

| Exposed Unexposed | Total

+
+

Cases | 32 7 | 39
Noncases | 28 53 | 81

+
+

Total | 60 60 | 120
| |
Risk | .5333333 .1166667 | .325
| |
| Point estimate | [95% Conf. Interval]
+
+

Risk difference | .4166667 | .2665567 .5667766
Risk ratio | 4.571429 | 2.191203 9.53721
Attr. frac. ex. | .78125 | .5436296 .8951475
Attr. frac. pop | .6410256 |
Odds ratio | 8.653061 | 3.441938 21.63764 (Cornfield)
+
+

chi2(1) = 23.74 Pr>chi2 = 0.0000

. generate dose2=(dose>0.16) if dose~=.

. tabulate dead dose2, chi2 column exact lrchi2 /* 2x2 table analysis */

+
+
| Key |
+
+
| frequency |
| column percentage |
+
+

| dose2
dead | 0 1 | Total

+
+

0 | 53 28 | 81
| 88.33 46.67 | 67.50

+
+

1 | 7 32 | 39
| 11.67 53.33 | 32.50

+
+

Total | 60 60 | 120
| 100.00 100.00 | 100.00

Pearson chi2(1) = 23.7417 Pr = 0.000
likelihood-ratio chi2(1) = 25.2010 Pr = 0.000
Fisher's exact = 0.000
1-sided Fisher's exact = 0.000

. logit dead dose2 /* logistic regression */

Iteration 0: log likelihood = -75.669723
Iteration 1: log likelihood = -63.522292
Iteration 2: log likelihood = -63.070143
Iteration 3: log likelihood = -63.06923
Iteration 4: log likelihood = -63.06923

Logistic regression Number of obs = 120
LR chi2(1) = 25.20
Prob > chi2 = 0.0000
Log likelihood = -63.06923 Pseudo R2 = 0.1665

+

--------

dead | Coef. Std. Err. z P>|z| [95% Con f. Interval]

+

--------

dose2 | 2.157913 .4782148 4.51 0.000 1.220629 3.095197
_cons | -2.024382 .4021506 -5.03 0.000 -2.812582 -1.236181

+

--------

. logistic dead dose2 /* same as logit with option or */

Logistic regression Number of obs = 120
LR chi2(1) = 25.20
Prob > chi2 = 0.0000
Log likelihood = -63.06923 Pseudo R2 = 0.1665

+

--------

dead | Odds Ratio Std. Err. z P>|z| [95% Con f. Interval]

+

--------

dose2 | 8.653061 4.138022 4.51 0.000 3.38932 22.09159
_cons | .1320755 .0531142 -5.03 0.000 .0600497 .2904914

+

--------

. display 1/(1+exp(2.024382))
.11666664

. logit dead dose

Iteration 0: log likelihood = -75.669723
Iteration 1: log likelihood = -64.201632
Iteration 2: log likelihood = -63.944944
Iteration 3: log likelihood = -63.944713
Iteration 4: log likelihood = -63.944713

Logistic regression Number of obs = 120
LR chi2(1) = 23.45
Prob > chi2 = 0.0000
Log likelihood = -63.944713 Pseudo R2 = 0.1549

+

--------

dead | Coef. Std. Err. z P>|z| [95% Con f. Interval]

+

--------

dose | 14.6369 3.332533 4.39 0.000 8.10526 21.16855
_cons | -3.569738 .7053432 -5.06 0.000 -4.952185 -2.187291

+

--------

. logistic dead dose /* same as logit with option or */

Logistic regression Number of obs = 120
LR chi2(1) = 23.45
Prob > chi2 = 0.0000
Log likelihood = -63.944713 Pseudo R2 = 0.1549

+

--------

dead | Odds Ratio Std. Err. z P>|z| [95% Con f. Interval]

+

--------

dose | 2273668 7577074 4.39 0.000 3311.843 1.56e+09
_cons | .0281632 .0198647 -5.06 0.000 .0070679 .1122204

+

--------

. display 1/(1+exp(3.569738))
.02739179

.
. * demonstration of analysis using grouped data
. use micegrp.dta, clear

. blogit dead total dose

Logistic regression for grouped data Number of obs = 120
LR chi2(1) = 23.45
Prob > chi2 = 0.0000
Log likelihood = -63.944713 Pseudo R2 = 0.1549

+

--------

_outcome | Coef. Std. Err. z P>|z| [95% Con f. Interval]

+

--------

dose | 14.6369 3.332533 4.39 0.000 8.10526 21.16855
_cons | -3.569738 .7053432 -5.06 0.000 -4.952185 -2.187291

+

--------

. predict linpred, xb

. predict phat, p

. predict deadhat, n

.
end of do-file

. exit, clear