library(finalfit)
ob =read.csv("~/Dropbox/_Conferences and Workshops/Cho Ray Workshop 2019/Datasets/obesity data.csv")
head(ob)
##   id gender height weight  bmi age  bmc  bmd   fat  lean pcfat
## 1  1      F    150     49 21.8  53 1312 0.88 17802 28600  37.3
## 2  2      M    165     52 19.1  65 1309 0.84  8381 40229  16.8
## 3  3      F    157     57 23.1  64 1230 0.84 19221 36057  34.0
## 4  4      F    156     53 21.8  56 1171 0.80 17472 33094  33.8
## 5  5      M    160     51 19.9  54 1681 0.98  7336 40621  14.8
## 6  6      F    153     47 20.1  52 1358 0.91 14904 30068  32.2

Simple description of data

xvars = c("age", "height", "weight", "pcfat")
yvar = c("gender")
ob %>% summary_factorlist(yvar, xvars,p=TRUE, add_dependent_label=TRUE) -> t1
knitr::kable(t1, row.names=FALSE)
Dependent: gender F M p
age Mean (SD) 48.6 (16.4) 43.7 (18.8) <0.001
height Mean (SD) 153.3 (5.6) 165.1 (6.7) <0.001
weight Mean (SD) 52.3 (7.7) 62.0 (9.6) <0.001
pcfat Mean (SD) 34.7 (5.2) 24.2 (5.8) <0.001
# Linear regression analysis
xvars = c("gender", "age", "bmi")
yvar = c("pcfat")
ob %>% finalfit(yvar, xvars) -> t2
## Dependent is not a factor and will be treated as a continuous variable
knitr::kable(t2, row.names=FALSE, align=c("l", "l", "r", "r", "r"))
Dependent: pcfat Mean (sd) Coefficient (univariable) Coefficient (multivariable)
gender F 34.7 (5.2) - -
M 24.2 (5.8) -10.52 (-11.18 to -9.85, p<0.001) -10.81 (-11.30 to -10.31, p<0.001)
age [13,88] 31.6 (7.2) 0.13 (0.11 to 0.15, p<0.001) 0.05 (0.03 to 0.06, p<0.001)
bmi [14.5,37.1] 31.6 (7.2) 1.04 (0.92 to 1.15, p<0.001) 1.09 (1.01 to 1.16, p<0.001)

Logistic regression analysis

os = read.csv("~/Dropbox/_Conferences and Workshops/Cho Ray Workshop 2019/Datasets/Hip fracture data.csv")
os$hipfx = as.factor(os$hipfx)
os$v4 = os$v4 / 0.15
head(os)
##   id     dov gender age      dob visit   v1   v2    v3       v4 wt bmi  ht
## 1  3 15/6/89   Male  73   8/6/16     1 0.98 0.88 1.079 9.720000 98  32 175
## 2  8 17/4/89 Female  67 11/12/21     1 0.85 0.85 0.966 8.833333 72  26 166
## 3  9 12/6/90   Male  68   8/1/22     1 0.87 0.84 1.013 9.960000 87  26 184
## 4 10  4/6/90 Female  62  15/5/28     1 0.62 0.71 0.839 8.093333 72  24 173
## 5 23  8/8/89   Male  61  22/9/28     1 0.87 0.60 0.811 7.626667 72  24 173
## 6 24  3/5/89 Female  76   1/8/13     1 0.76 0.58 0.743 6.533333 67  28 156
##   v5   v6 v7 v8 v9 hipfx timehip
## 1 NA 39.9  1  0  0     0    0.55
## 2 18 31.0  0  0  0     0   19.68
## 3 36 28.6  0  0  0     0    5.05
## 4 NA 28.2  1  0  0     0   18.55
## 5 44 28.9  1  0  0     0   19.37
## 6 15 33.3  0  0  0     0   12.30
xvars = c("gender", "age", "bmi", "v4")
yvar = c("hipfx")
os %>% finalfit(yvar, xvars) -> t3
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
knitr::kable(t3, row.names=FALSE, align=c("l", "l", "r", "r", "r"))
Dependent: hipfx 0 1 OR (univariable) OR (multivariable)
gender Female 1512 (58.2) 142 (75.1) - -
Male 1087 (41.8) 47 (24.9) 0.46 (0.32-0.64, p<0.001) 0.92 (0.61-1.36, p=0.676)
age Mean (SD) 68.9 (6.3) 75.3 (7.7) 1.13 (1.10-1.15, p<0.001) 1.10 (1.08-1.12, p<0.001)
bmi Mean (SD) 26.9 (4.7) 23.9 (3.9) 0.84 (0.81-0.88, p<0.001) 0.91 (0.87-0.96, p<0.001)
v4 Mean (SD) 7.8 (1.5) 6.6 (1.5) 0.55 (0.49-0.62, p<0.001) 0.67 (0.58-0.77, p<0.001)
# Plotting odds ratio
os %>% 
  or_plot(yvar, xvars)
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Warning: Removed 1 rows containing missing values (geom_errorbarh).