Introduction

Markdown

This is an R Markdown html notebookdocument. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.

require(rms)

Data

Setup

getHdata(titanic3)  # Get the dataset from the VU DataSets page
mu <- markupSpecs$html   # markupSpecs is in Hmisc
subtext <- mu$subtext
code    <- mu$code

Data Dictionary

html(contents(titanic3), maxlevels=10, levelType='table')

Data frame:titanic3

1309 observations and 14 variables, maximum # NAs:1188  
NameLabelsUnitsLevelsStorageNAs
pclass 3integer 0
survivedSurvivedinteger 0
nameNamecharacter 0
sex 2integer 0
ageAgeYeardouble 263
sibspNumber of Siblings/Spouses Aboardinteger 0
parchNumber of Parents/Children Aboardinteger 0
ticketTicket Numbercharacter 0
farePassenger FareBritish Pound (\243)double 1
cabin187integer 0
embarked 3integer 2
boat 28integer 0
bodyBody Identification Numberinteger1188
home.destHome/Destinationcharacter 0

VariableLevels
pclass1st
2nd
3rd
sexfemale
male
cabin
A10
A11
A14
A16
A18
A19
A20
A21
A23
...
embarkedCherbourg
Queenstown
Southampton
boat
1
10
11
12
13
13 15
13 15 B
14
15
...

Descriptive Statistics
for the titanic3 dataset

# Set graphics type so that Hmisc and rms packages use plotly
# Chunk header height=150 is in pixels
# For certain print methods set to use html
options(grType='plotly', prType='html')
s <- summaryM(age + pclass ~ sex, data=titanic3)
s


Descriptive Statistics  (N=1309)

+------------+----+-----------+-----------+
|            |N   |female     |male       |
|            |    |(N=466)    |(N=843)    |
+------------+----+-----------+-----------+
|Age [Year]  |1046|  19/27/38 |  21/28/39 |
+------------+----+-----------+-----------+
|pclass : 1st|1309|0.31  (144)|0.21  (179)|
+------------+----+-----------+-----------+
|    2nd     |    |0.23  (106)|0.20  (171)|
+------------+----+-----------+-----------+
|    3rd     |    |0.46  (216)|0.58  (493)|
+------------+----+-----------+-----------+
plot(s)
Ignoring 2 observations
$Categorical


$Continuous
d <- describe(titanic3)
plot(d)
$Categorical
Numeric color variables cannot (yet) be mapped to text.
Feel free to make a feature request 
https://github.com/plotly/plotly.jsNumeric color variables cannot (yet) be mapped to text.
Feel free to make a feature request 
https://github.com/plotly/plotly.js


$Continuous

NA

The following doesn’t work because it overlays two different legends

# Try combining two plots into one
p <- plot(d)
plotly::subplot(p[[1]], p[[2]],
                nrows=2, heights=c(.3, .7), which_layout=1)

Logistic Regression Model

dd <- datadist(titanic3); options(datadist='dd')
f <- lrm(survived ~ rcs(sqrt(age),5) * sex, data=titanic3)
f
latex(f)

\[{\rm Prob}\{{\rm survived}=1\} = \frac{1}{1+\exp(-X\beta)}, {\rm \ \ where} \\ \] \begin{eqnarray*} X\hat{\beta}= & & \\ & & 0.6476552 \\ & & + 0.03027852 {\rm age}+0.01114664 ({\rm age}-2.236068)_{+}^{3}-0.08633279({\rm age}-4.582576)_{+}^{3} \\ & & +0.1857246 ({\rm age}-5.291503)_{+}^{3}-0.1516536 ({\rm age}-6.082763)_{+}^{3} \\ & & +0.04111512 ({\rm age}-7.549834)_{+}^{3} \\ & & +1.345826[{\rm male}] \\ & & +[{\rm male}][-1.036388 {\rm age}+0.06620937 ({\rm age}-2.236068)_{+}^{3} \\ & & -0.5363316 ({\rm age}-4.582576)_{+}^{3}+0.6312123 ({\rm age}-5.291503)_{+}^{3} \\ & & -0.1266968 ({\rm age}-6.082763)_{+}^{3}-0.03439327({\rm age}-7.549834)_{+}^{3} ] \\ \end{eqnarray*} and \([c]=1\) if subject is in group \(c\), 0 otherwise; \((x)_{+}=x\) if \(x > 0\), 0 otherwise
\({\rm age}\) is pre--transformed as \({\rm \sqrt{age}}\).

a <- anova(f)
a
plot(a)
s <- summary(f, age=c(2, 21))
plot(s, log=TRUE)

print(s, dec=2)
ggplot(Predict(f, age, sex), height=500, width=650)  # uses ggplotly()

plotp(Predict(f, age, sex))                          # uses plotly directly
minimal value for n is 3, returning requested palette with 3 different levels
minimal value for n is 3, returning requested palette with 3 different levels

plot(nomogram(f, fun=plogis, funlabel='Prob(survive)'))

Survival Plots for pbc Dataset

getHdata(pbc)
pbc <- upData(pbc, 
              fu.yrs = fu.days / 365.25,
              units = c(fu.yrs = 'year'))
Input object size:   76592 bytes;    19 variables    418 observations
Added variable      fu.yrs
New object size:    80712 bytes;    20 variables    418 observations
f <- npsurv(Surv(fu.yrs, status) ~ spiders, data=pbc)
survplotp(f, time.inc=1, times=c(5, 10))

Computing Environment

 R version 3.3.3 (2017-03-06)
 Platform: x86_64-pc-linux-gnu (64-bit)
 Running under: Ubuntu 16.04.2 LTS
 
 attached base packages:
 [1] stats     graphics  grDevices utils     datasets 
 [6] methods   base     
 
 other attached packages:
 [1] rms_5.1-1       SparseM_1.76    Hmisc_4.0-3    
 [4] ggplot2_2.2.1   Formula_1.2-1   survival_2.41-3
 [7] lattice_0.20-35
 
To cite R in publication use:

R Core Team (2017). R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria. https://www.R-project.org/.

LS0tCnRpdGxlOiAiRXhhbXBsZXMgZm9yIHJtcyBQYWNrYWdlIgphdXRob3I6ICJGRSBIYXJyZWxsIgpkYXRlOiAnYHIgU3lzLkRhdGUoKWAnCm91dHB1dDoKICBodG1sX25vdGVib29rOgogICAgdG9jOiB5ZXMKICAgIHRvY19mbG9hdDoKICAgICAgY29sbGFwc2VkOiB5ZXMKICBodG1sX2RvY3VtZW50OgogICAgdG9jOiB5ZXMKLS0tCiMgSW50cm9kdWN0aW9uCiMjIE1hcmtkb3duClRoaXMgaXMgYW4gUiBNYXJrZG93biBodG1sIG5vdGVib29rZG9jdW1lbnQuIE1hcmtkb3duIGlzIGEgc2ltcGxlIGZvcm1hdHRpbmcgc3ludGF4IGZvciBhdXRob3JpbmcgSFRNTCwgUERGLCBhbmQgTVMgV29yZCBkb2N1bWVudHMuIEZvciBtb3JlIGRldGFpbHMgb24gdXNpbmcgUiBNYXJrZG93biBzZWUgPGh0dHA6Ly9ybWFya2Rvd24ucnN0dWRpby5jb20+LgoKYGBge3IsIHJlc3VsdHM9J2hpZGUnfQpyZXF1aXJlKHJtcykKYGBgCgojIERhdGEgey50YWJzZXR9CiMjIFNldHVwCmBgYHtyIHQzfQpnZXRIZGF0YSh0aXRhbmljMykgICMgR2V0IHRoZSBkYXRhc2V0IGZyb20gdGhlIFZVIERhdGFTZXRzIHBhZ2UKbXUgPC0gbWFya3VwU3BlY3MkaHRtbCAgICMgbWFya3VwU3BlY3MgaXMgaW4gSG1pc2MKc3VidGV4dCA8LSBtdSRzdWJ0ZXh0CmNvZGUgICAgPC0gbXUkY29kZQpgYGAKIyMgRGF0YSBEaWN0aW9uYXJ5CmBgYHtyIGRkaWN0fQpodG1sKGNvbnRlbnRzKHRpdGFuaWMzKSwgbWF4bGV2ZWxzPTEwLCBsZXZlbFR5cGU9J3RhYmxlJykKYGBgCgojIyBEZXNjcmlwdGl2ZSBTdGF0aXN0aWNzYHIgc3VidGV4dCgnZm9yIHRoZScsIGNvZGUoJ3RpdGFuaWMzJyksICdkYXRhc2V0JylgCmBgYHtyIHQzZCwgaGVpZ2h0PTE1MH0KIyBTZXQgZ3JhcGhpY3MgdHlwZSBzbyB0aGF0IEhtaXNjIGFuZCBybXMgcGFja2FnZXMgdXNlIHBsb3RseQojIENodW5rIGhlYWRlciBoZWlnaHQ9MTUwIGlzIGluIHBpeGVscwojIEZvciBjZXJ0YWluIHByaW50IG1ldGhvZHMgc2V0IHRvIHVzZSBodG1sCm9wdGlvbnMoZ3JUeXBlPSdwbG90bHknLCBwclR5cGU9J2h0bWwnKQpzIDwtIHN1bW1hcnlNKGFnZSArIHBjbGFzcyB+IHNleCwgZGF0YT10aXRhbmljMykKcwpwbG90KHMpCmQgPC0gZGVzY3JpYmUodGl0YW5pYzMpCnBsb3QoZCkKYGBgCgpUaGUgZm9sbG93aW5nIGRvZXNuJ3Qgd29yayBiZWNhdXNlIGl0IG92ZXJsYXlzIHR3byBkaWZmZXJlbnQgbGVnZW5kcwpgYGB7ciBzdWIsaGVpZ2h0PTYwMCxldmFsPUZBTFNFfQojIFRyeSBjb21iaW5pbmcgdHdvIHBsb3RzIGludG8gb25lCnAgPC0gcGxvdChkKQpwbG90bHk6OnN1YnBsb3QocFtbMV1dLCBwW1syXV0sCiAgICAgICAgICAgICAgICBucm93cz0yLCBoZWlnaHRzPWMoLjMsIC43KSwgd2hpY2hfbGF5b3V0PTEpCmBgYAoKIyBMb2dpc3RpYyBSZWdyZXNzaW9uIE1vZGVsCmBgYHtyIGxybXR9CmRkIDwtIGRhdGFkaXN0KHRpdGFuaWMzKTsgb3B0aW9ucyhkYXRhZGlzdD0nZGQnKQpmIDwtIGxybShzdXJ2aXZlZCB+IHJjcyhzcXJ0KGFnZSksNSkgKiBzZXgsIGRhdGE9dGl0YW5pYzMpCmYKbGF0ZXgoZikKYSA8LSBhbm92YShmKQphCnBsb3QoYSkKYGBgCgpgYGB7ciBzdW1tYXJ5fQpzIDwtIHN1bW1hcnkoZiwgYWdlPWMoMiwgMjEpKQpwbG90KHMsIGxvZz1UUlVFKQpwcmludChzLCBkZWM9MikKYGBgCgpgYGB7ciBnZ3AsZmlnLmhlaWdodD01LGZpZy53aWR0aD02fQpnZ3Bsb3QoUHJlZGljdChmLCBhZ2UsIHNleCksIGhlaWdodD01MDAsIHdpZHRoPTY1MCkgICMgdXNlcyBnZ3Bsb3RseSgpCnBsb3RwKFByZWRpY3QoZiwgYWdlLCBzZXgpKSAgICAgICAgICAgICAgICAgICAgICAgICAgIyB1c2VzIHBsb3RseSBkaXJlY3RseQpwbG90KG5vbW9ncmFtKGYsIGZ1bj1wbG9naXMsIGZ1bmxhYmVsPSdQcm9iKHN1cnZpdmUpJykpCmBgYAoKIyBTdXJ2aXZhbCBQbG90cyBmb3IgYHIgbXUkY29kZSgncGJjJylgIERhdGFzZXQKYGBge3IgcGJjLGZpZy5oZWlnaHQ9NixmaWcud2lkdGg9N30KZ2V0SGRhdGEocGJjKQpwYmMgPC0gdXBEYXRhKHBiYywgCiAgICAgICAgICAgICAgZnUueXJzID0gZnUuZGF5cyAvIDM2NS4yNSwKICAgICAgICAgICAgICB1bml0cyA9IGMoZnUueXJzID0gJ3llYXInKSkKZiA8LSBucHN1cnYoU3VydihmdS55cnMsIHN0YXR1cykgfiBzcGlkZXJzLCBkYXRhPXBiYykKc3VydnBsb3RwKGYsIHRpbWUuaW5jPTEsIHRpbWVzPWMoNSwgMTApKQpgYGAKCiMgQ29tcHV0aW5nIEVudmlyb25tZW50CmByIG11JHNlc3Npb24oKWAKCgo8IS0tLQpUbyB1cGRhdGUgaHRtbCBub3RlYm9vayBvbiB0aGUgc2VydmVyOiBjZGF0YXIgZXhhbXBsZXMubmIuaHRtbCBybXMKLS0+Cg==