library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.2 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ ggplot2 3.4.2 ✔ tibble 3.2.1
## ✔ lubridate 1.9.2 ✔ tidyr 1.3.0
## ✔ purrr 1.0.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(dplyr)
library(nnet)
library(readstata13)
duc = read.dta13("D:\\OneDrive - UMP\\Desktop\\ducyhdp.dta")
## Warning in read.dta13("D:\\OneDrive - UMP\\Desktop\\ducyhdp.dta"):
## Factor codes of type double or float detected in variables
##
## a1, a4, b22, plnhanthuc, plnhanthuc1
##
## No labels have been assigned.
## Set option 'nonint.factors = TRUE' to assign labels anyway.
## Warning in read.dta13("D:\\OneDrive - UMP\\Desktop\\ducyhdp.dta"):
## Missing factor labels for variables
##
## c3_3
##
## No labels have been assigned.
## Set option 'generate.factors=TRUE' to generate labels.
attach(duc)
names(duc)
## [1] "a1" "a2" "a3" "a4" "a5"
## [6] "a6" "b1" "b2" "b3" "b4"
## [11] "b5" "b6" "b7" "b8" "b9"
## [16] "b10" "b11" "b12" "b14" "b15"
## [21] "b16" "b17" "b18" "b19" "b20"
## [26] "b21" "b22" "b23" "b24" "b25"
## [31] "b26" "b27" "c1_1" "c1_2" "c1_3"
## [36] "c1_4" "c1_5" "c1_6" "c2" "c3_1"
## [41] "c3_2" "c3_3" "c3_4" "B1dung" "B2dung"
## [46] "B3dung" "B6dung" "B17dung" "B18dung" "B20dung"
## [51] "B13" "nhanthuc" "plnhanthuc" "plnhanthuc1"
# Câu lệnh trong stata - glm phu_thuoc i.doc_lap, family(poisson) robust eform
# glm plnhanthuc1 i.a3, family(poisson) robust eform
##Iteration 0: log pseudolikelihood = -236.30018
##Iteration 1: log pseudolikelihood = -227.73641
##Iteration 2: log pseudolikelihood = -227.67647
##Iteration 3: log pseudolikelihood = -227.67644
##Iteration 4: log pseudolikelihood = -227.67644
##Generalized linear models No. of obs = 414
##Optimization : ML Residual df = 412
## Scale parameter = 1
##Deviance = 253.3528876 (1/df) Deviance = .6149342
##Pearson = 272.9769185 (1/df) Pearson = .6625653
##
##Variance function: V(u) = u [Poisson]
##Link function : g(u) = ln(u) [Log]
## AIC = 1.109548
##Log pseudolikelihood = -227.6764438 BIC = -2229.304
##------------------------------------------------------------------------------
## | Robust
## plnhanthuc1 | IRR Std. Err. z P>|z| [95% Conf. Interval]
##-------------+----------------------------------------------------------------
## a3 | .5207694 .0392859 -8.65 0.000 .4491924 .603752
## _cons | .8950968 .1402836 -0.71 0.479 .6583636 1.216954
##------------------------------------------------------------------------------
# glm plnhanthuc1 i.a3, family(poisson) robust eform
##Iteration 0: log pseudolikelihood = -215.88496
##Iteration 1: log pseudolikelihood = -207.77774
##Iteration 2: log pseudolikelihood = -207.26174
##Iteration 3: log pseudolikelihood = -207.13962
##Iteration 4: log pseudolikelihood = -207.11063
##Iteration 5: log pseudolikelihood = -207.10448
##Iteration 6: log pseudolikelihood = -207.1032
##Iteration 7: log pseudolikelihood = -207.10291
##Iteration 8: log pseudolikelihood = -207.10284
##Iteration 9: log pseudolikelihood = -207.10282
##
##Generalized linear models No. of obs = 414
##Optimization : ML Residual df = 410
## Scale parameter = 1
##Deviance = 212.2056463 (1/df) Deviance = .5175747
##Pearson = 267.0153185 (1/df) Pearson = .6512569
##
##Variance function: V(u) = u [Poisson]
##Link function : g(u) = ln(u) [Log]
##
## AIC = 1.01982
## pseudolikelihood = -207.1028231 BIC = -2258.399
##
##------------------------------------------------------------------------------
## | Robust
## plnhanthuc1 | IRR Std. Err. z P>|z| [95% Conf. Interval]
##-------------+----------------------------------------------------------------
## a3 |
## kha | 1.430492 .2585881 1.98 0.048 1.003726 2.038713
## trung binh | .1194788 .0616065 -4.12 0.000 .0434899 .3282414
## yeu | 2.33e-07 5.00e-08 -71.18 0.000 1.53e-07 3.55e-07
## |
## _cons | .298931 .0465394 -7.76 0.000 .2203184 .4055938
##------------------------------------------------------------------------------
Cau lenh trong R
library(biostat3) # eform
## Loading required package: survival
## Loading required package: MASS
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
##
## Attaching package: 'biostat3'
## The following object is masked from 'package:survival':
##
## colon
## The following object is masked from 'package:lubridate':
##
## year
library(car)
## Loading required package: carData
##
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
##
## recode
## The following object is masked from 'package:purrr':
##
## some
m1 = glm(plnhanthuc1 ~ a3)
eform(m1)
## exp(beta) 2.5 % 97.5 %
## (Intercept) 1.3484679 1.2471933 1.4579662
## a3kha 1.1373533 1.0300815 1.2557963
## a3trung binh 0.7685461 0.6907990 0.8550433
## a3yeu 0.7415823 0.6462143 0.8510248
lincom(m1,c("a3kha", "a3trung binh","a3yeu"),eform=TRUE)
## Estimate 2.5 % 97.5 % Chisq Pr(>Chisq)
## a3kha 1.137353 1.030081 1.255796 6.483811 0.01088613
## a3trung binh 0.7685461 0.690799 0.8550433 23.4053 1.312166e-06
## a3yeu 0.7415823 0.6462143 0.8510248 18.12025 2.07384e-05
exp(confint(m1))
## Waiting for profiling to be done...
## 2.5 % 97.5 %
## (Intercept) 1.2471933 1.4579662
## a3kha 1.0300815 1.2557963
## a3trung binh 0.6907990 0.8550433
## a3yeu 0.6462143 0.8510248
exp(coef(m1))
## (Intercept) a3kha a3trung binh a3yeu
## 1.3484679 1.1373533 0.7685461 0.7415823
het