转自医学方
2019-07-4 Alexander
流行病学或者医学论文中，对研究对象基本情况的描述通常以表格的形式进行，并且放在结果部分的开头，即Table 1，主要内容是研究对象一般情况和研究变量或协变量的分组展示。
前几天文章修回过程中，花了两天时间分析数据，修改文章，其中有近1天的时间都在手动录入数据（从R studio里把分析结果整理到Excel或者word），这样除了花费时间外，还非常容易出错。之前一直想找时间通过R markdown把制作表格的过程程序化，可是效果并不理想。
这次痛定思痛，先从table 1开始，发现了几个不错的方法。其中一种个人觉得可读性和可编辑性都比较强，于是学习了一下，作为一个非常实用的工具分享给大家。
这里主要参考一篇博客Fast-track publishing using knitr: table mania，对细节进行了加工和注释。
1 数据的准备
数据主要来自于boot包的melanoma。加载后，看下数据的基本结构。
接下来对数据进行简单的整理，为后续分析做准备；
将分类变量定义为因子型并设置标签（这里建议设置一个新的变量，仅用于table 1的制作，不影响后续的分析）；
2 安装和加载R包 Gmisc
后面两个包是加载“Gmisc”时要求加载的。
3 自定义函数、制作表格
根据已有函数自定义函数，并制作表格。定义一个函数，输入数据集的变量并得到该变量的统计结果：
函数定义完成后，建立一个空的列表，以储存每个变量的分析结果，并进行分析，将结果储存在列表中：
将所有结果merge到一个矩阵中，并建立rgroup（table1第一列的变量名） 和 n.rgroup（table 1第一列每个变量的行数）：
结果如下：
当然，有些情况下，需要多加一个分组标题栏（column spanner），该怎么加呢？
如下：
结果如下：
4 导出结果
在R studio viewer窗口点击白色按钮，即可在浏览器中打开，然后复制粘贴到word可以进一步加工修饰。
是不是很刺激呢。
应该还有其他的导出方法，不过这个已经很方便了。
拓展功能选
⒈ 二分类变量只显示一个（比如男性和女性）。只要在getDescriptionStatsBy的"show_all_values"参数设置为FALSE即可；
⒉ 显示缺失值。getDescriptionStatsBy的"useNA"参数设置为"ifany"，表示如果有缺失值就显示缺失值情况；如设置为“no”，表示始终不显示缺失值情况；“always”则表示无论是否有缺失值都显示缺失值情况；
⒊ Total一列是可以去掉的，getDescriptionStatsBy的"add_total_col"参数设置为FALSE即可。
不足之处
⒈ 差异性检验是采用非参的方法，虽然没有错，但是一般符合参数检验条件的数据还是要使用参数检验的方法，这里可以自行检验后再修改P-value；
⒉ Mean (SD)的展示形式有个括号感觉有点别扭，还不知道怎么去掉，有方法的小伙伴欢迎分享交流。
另外有一些其他的制作table 1的R包，比如table 1（R包的名字）包，tableone包，还有其他生成表格的R包（plyr等），个人浏览下来感觉这个最容易理解和掌握，其他包的功能有兴趣的可以再自行挖掘对比。
原文链接： https://mp.weixin.qq.com/s?src=11×tamp=1562230826&ver=1707&signature=Og8lYPNfFi99QvnQb8OAtkTIo75N9G0JHqvpXxLS5aRpqRcnlgtxXJAMtfgxB8kAK8vinKSxdO6A1qxNy-4k8AyE9wUMYKSarBLDydWO-vazmCNPJIAa5GfaBiFIghaO&new=1
数据分析之美：决策树R语言实现R语言实现决策树
1.准备数据
[plain] view plain copy
>install.packages("tree")
>library(tree)
>library(ISLR)
>attach(Carseats)
>High=ifelse(Sales<=8,"No","Yes") //set high values by sales data to calssify
>Carseats=data.frame(Carseats,High) //include the high data into the data source
>fix(Carseats)
2.生成决策树
[plain] view plain copy
>tree.carseats=tree(High~.-Sales,Carseats)
>summary(tree.carseats)
[plain] view plain copy
//output training error is 9%
Classification tree:
tree(formula = High ~ . - Sales, data = Carseats)
Variables actually used in tree construction:
[1] "ShelveLoc" "Price" "Income" "CompPrice" "Population"
[6] "Advertising" "Age" "US"
Number of terminal nodes: 27
Residual mean deviance: 0.4575 = 170.7 / 373
Misclassification error rate: 0.09 = 36 / 400
3. 显示决策树
[plain] view plain copy
>plot(tree . carseats )
>text(tree .carseats ,pretty =0)
4.Test Error
[plain] view plain copy
//prepare train data and test data
//We begin by using the sample() function to split the set of observations sample() into two halves, by selecting a random subset of 200 observations out of the original 400 observations.
>set . seed (1)
>train=sample(1:nrow(Carseats),200)
>Carseats.test=Carseats[-train,]
>High.test=High[-train]
//get the tree model with train data
>tree. carseats =tree (High~.-Sales , Carseats , subset =train )
//get the test error with tree model, train data and predict method
//predict is a generic function for predictions from the results of various model fitting functions.
>tree.pred = predict ( tree.carseats , Carseats .test ,type =" class ")
>table ( tree.pred ,High. test)
High. test
tree. pred No Yes
No 86 27
Yes 30 57
>(86+57) /200
[1] 0.715
5.决策树剪枝
[plain] view plain copy
/**
Next, we consider whether pruning the tree might lead to improved results. The function cv.tree() performs cross-validation in order to cv.tree() determine the optimal level of tree complexitycost complexity pruning is used in order to select a sequence of trees for consideration.
For regression trees, only the default, deviance, is accepted. For classification trees, the default is deviance and the alternative is misclass (number of misclassifications or total loss).
We use the argument FUN=prune.misclass in order to indicate that we want the classification error rate to guide the cross-validation and pruning process, rather than the default for the cv.tree() function, which is deviance.
If the tree is regression tree,
>plot(cv. boston$size ,cv. boston$dev ,type=’b ’)
*/
>set . seed (3)
>cv. carseats =cv. tree(tree .carseats ,FUN = prune . misclass ,K=10)
//The cv.tree() function reports the number of terminal nodes of each tree considered (size) as well as the corresponding error rate(dev) and the value of the cost-complexity parameter used (k, which corresponds to α.
>names (cv. carseats )
[1] " size" "dev " "k" " method "
>cv. carseats
$size //the number of terminal nodes of each tree considered
[1] 19 17 14 13 9 7 3 2 1
$dev //the corresponding error rate
[1] 55 55 53 52 50 56 69 65 80
$k // the value of the cost-complexity parameter used
[1] -Inf 0.0000000 0.6666667 1.0000000 1.7500000
2.0000000 4.2500000
[8] 5.0000000 23.0000000
$method //miscalss for classification tree
[1] " misclass "
attr (," class ")
[1] " prune " "tree. sequence "
[plain] view plain copy
//plot the error rate with tree node size to see whcih node size is best
>plot(cv. carseats$size ,cv. carseats$dev ,type=’b ’)
/**
Note that, despite the name, dev corresponds to the cross-validation error rate in this instance. The tree with 9 terminal nodes results in the lowest cross-validation error rate, with 50 cross-validation errors. We plot the error rate as a function of both size and k.
*/
>prune . carseats = prune . misclass ( tree. carseats , best =9)
>plot( prune . carseats )
>text( prune .carseats , pretty =0)
//get test error again to see whether the this pruned tree perform on the test data set
>tree.pred = predict ( prune . carseats , Carseats .test , type =" class ")
>table ( tree.pred ,High. test)
High. test
tree. pred No Yes
No 94 24
Yes 22 60
>(94+60) /200
[1] 0.77
R语言数据对象与运算R语言数据对象与运算 笔记整理2.1 数据对象及类型R语言创建和控制的实体被称为对象（object）ls()命令来查看当前系统里的数据对象R对象的名称必须以一个英文字母打头，并由一串大小写字母、数字或钟点组成注意：R区分大小写不要用R的内置函数名称作为数据对象的名称，如c、length等2.2 数据对象类型R语言的对象包括数值型（numeric）：实数， 可写成整数（integers）、小数（decimal fractions）、科学记数（scientific notation）逻辑型（logical）：T(true)或F（FALSE）字符型（character）：夹在" "或之间复数型（complex）：形如a+bi原味型（raw）：以二进制形式保存数据缺省型（missing value）：有些统计资料是不完整的，当一个元素或值在统计的时候是“不可得到（not available）”或“缺失值（missing value）”的时候，相关位置可能会被保留并赋予一个特定的NA（not available）值，任何NA的运算结果都是NA。辨别和转换数据对象类型的函数：辨别 转换character is.character() as,character()complexdoubleintegerlogicalNAnumeric2.3 数据对象构造R语言里的数据对象主要有六种构造：向量（vector）、矩阵（matrix）、数组（array）、列表（list）、数据框（data frames）、因子（factor）2.3.1 向量（vector）是由有相同基本类型元素组成的序列，相当于一维数组 5个数值组成的向量x，这是一个用函数c()完成的赋值语句，这里c()可以有任意多个参数，而它输出的值则是一个把这些参数首尾相连形成的一个向量R的赋值符号除了“<-”外，还有"->""="例如：>c(1,3,5,7,9) ->y >y [1] 2 5 8 3 >z = c(1,3,5,7,9) >z [1] 1 3 5 7 9 assign（）函数对向量进行赋值 length（）：可返回向量的长度 mode（）可返回向量的数据类型 正则序列 用 “：”符号，可产生有规律的正则序列（： 的运算级别最高） 函数seq（）产生有规律的各种序列seq（from,to ,by） from 给序列的起始值，to表示序列的终止值，by表示步长（by 省略时，表示步长值为1）>seq(1,10,2) [1] 1 3 5 7 9 >seq(1,10) [1] 1 2 3 4 5 6 7 8 9 10 有时关注的是数列的长度，利用句法：seq（下界，by=，length=）>seq(1,by=2,length=10) [1] 1 3 5 7 9 11 13 15 17 19 rep(x,times，……)x表示要重复的对象，times表示重复的次数>rep(c(1,3),4) [1] 1 3 1 3 1 3 1 3 >rep(c(1,3),each=4) [1] 1 1 1 1 3 3 3 3 对每个元素进行重复R中的内置函数：mean（）来示向量的均值median（）求是位数var（）求方差sd（）求标准差sort（）对向量排序rev（）将向量按原方向的反方向排列rank（）给求出向量的秩prod（）求向量连乘积append（）为向量添加元素对向量运算常见函数表 函数 用途sum（） 求和max（） 求最大值min（） 求最小值range（） 求极差（全矩）mean（） 求均值median 求中位数var（） 求方差sd（） 求标准差sort（） 排序rev（） 反排序rank（） 求秩append（） 添加replace（） 替换match（） 匹配pmatch（） 部分匹配all（） 判断所有any（） 判断部分prod（） 积 2.3.2 矩阵矩阵（matrix)是将数据用行和列排列的长方形表格，它是二维的数组，其单元必须是相同的数据类型，通常用列来表示不同的变量，用行表示各个对象。其句法是：matrix（data=NA,ncol=1,byrow-=FALSE,dimnames=NULL）data是必须的，其它几个选择参数。nrow表示矩阵的行数ncol表示矩阵的列数byrow默认为FALSE,表示矩阵按列排列，如设置为T,表示按行排列；dimnames可更改矩阵行列名字diag()函数生成对角矩阵diag（）这个函数比较特别，当数据是向量时则生成对角矩阵，但当数据是矩阵时，则返回对角元素也可用函数diag（）生成单位矩阵 当我们生成了某个矩阵后，若要访问矩阵的某个元素或某行（列），可以利用形如A[i,j]的形式得到相应的索引矩阵矩阵可进行相应的加减乘除运算，但运算过程中要注意行数和列数的限制条件R里A*B并不是表示矩阵相乘，只表示矩阵对应的元素相乘矩阵相乘应用A%*%Bdim（）返回矩阵的行数和列数nrow（）返回矩阵的行数ncol（）返回矩阵的列数solve（）返回矩阵的逆矩阵对矩阵运算的常见函数 函数 用途as.matrix() 把非矩阵的转换成矩阵is.matrix() 辨别是否矩阵diag() 返回对角元素或生成对角矩阵eigen() 求特征值和特征向量solve() 求逆矩阵chol() Choleski分解svd() 奇异值分解qr() QR分解det() 求行列式dim() 返回行列数t() 矩阵转置apply() 对矩阵应用函数R语言还提供了专门针对矩阵的行或列计算的函数如 colSUms（）对矩阵各列求和colMeans（）求矩阵各列的均值类似的有 rowSums（）rowMeans（）更一般的方法：apply（）函数来对各行各列进行运算句法是：apply（X,MARGIN,FUN,……）X表示要处理的数据MARGIN表示函数作用的范围取1表示对行运用函数取2表示对列运用函数FUN表示要运用的函数rbind（）、cbind（）将两个或两个以上的矩阵合并起来rbind（）表示按行合并，cbind（）则表示按列合并2.3.3 数组数组（array）可以看作是带有多个下标的类型相同的元素的集合。数组的生成函数是array（），其句法是array（data=NA,dim=length(data),dimnames-NULL）data表示数据，可以为空dim 表示维数dimnames可以更改数组难度的名称2.3.4 列表向量、矩阵和的单元必须是同一类型的数据，若一个数据对象需要含有不同的数据类型，可采用列表（list）这种数据对象的形式。列表是一个对象的有序集合构成的对象，列表中包含的对象又称为它的分量（components），分量可以是不同的模式或（和）类型语法式为：list （变量1=分量1，变量2=分量2,……）若要访问列表的某一成分，可以用LST[[1]],LST[[2]]的形式访问因分量可以被命名，故可以在列表名称后加$符号，再写上成分名称来访问列表分量函数length（）、mode（）、names（）可以分别返回列表的长度（分量的数目）、数据类型、列表里成分的名字2.3.5 数据框数据框（data frame）是一种矩阵形式的数据，但数据框中各列可以是不同类型的数据。数据框每列是一个变量，每行是一个观测 。对可能列入数据框中的列表有如下的一些限制：1．分量必须是向量（数值，字符，逻辑），因子，数值矩阵，列表或者其他数据框。2．矩阵，列表和数据框为新的数据框提供了尽可能多的变量，因为它们各自拥有列、元素或者变量。3．数值向量、逻辑值、因子保持原有格式，而字符向量会被强制转换成因子并且它的水平就是向量中出现的独立值。4．在数据框中以变量形式出现的向量结构必须长度一致，矩阵结构必须有一样的行数。R中用函数data.frame()生成数据框，其句法是：data.frame（data1,data2,……）数据框的列名默认为变量名，也可对列名进行重新命名也可以对数据框的行名进行修改2.3.6 因子和有序因子分类型数据经常要把数据分成不同的水平或因子（factor）生成因子的命令是factor（），其句法是：factor（data,levels,labels,……）其中data表示数据levels是因子水平向量labels是因子的标签向量levels，labels是备选项，可以不选若上面的每个因子并不表示因子的大小，要表达因子之间有大小顺序（考虑因子之间的顺序），则可以用 ordered（）函数产生2.4 数据的录入及编辑c函数：c函数是把各个值联成一个向量或列表，可以形成数值型向量、字符型向量或其它类型向量 scan函数：功能类似于c函数，实际上是一种键盘输入数据函数。当输入scan（），然后按回车键，这时将等待输入数据，数据之间只要空格分开即可（c函数要用逗号分开）。输入完数据，再按回车键，这时数据录入完毕。scan函数还可以读入外部文本文件，若现有一个文本文件，data.txt，读入这个文件的命令是：>x=scan(file="dat.txt")若原文件的数据之间有逗号等分隔符，用scan读入应该去掉这些分隔符，其命令是：>x=scan(file="dat.txt"，sep=",") 编辑数据data.entry命令xx原先未被定义，现在赋予其一个空值，这时会出现一个电子表格界面，等待输入数据：>data.entry(xx=c(NA)) 当电子表格关闭后，数据会自动保存edit命令用来编辑函数，也可用来编辑数据，但不会自动保存fix函数与edit类似，但它可以自动保存从外部文件读入数据从文本文件读取：>s1=read.table("student.txt") >s1 V1V2V3 1class sexscore 2 1 女80 3 1 男85 4 2 男92 5 2 女76 6 3 女61 7 3 女95 8 3 男83 读入表格数据的命令是：read.table 忽略掉标签而直接使用默认的行标签>s2=read.table("student.txt",header=T)>s2 class sexscore 1 1 女80 2 1 男85 3 2 男92 4 2 女76 5 3 女61 6 3 女95 7 3 男83 从网络读入数据url可以从网页上读入正确格式的数据，要借助read.table函数> address=http://www.the-data-mine.com/bin/view/Misc/WebHome/sample.txt>read.table(file=url(address)) 读入其他格式的数据库要读入其他格式的数据库，必须先安装"foreign"模块，它不属于R的8个内置模块，需在使用前安装。 library(foreign) SAS：R只能诗篇SAS Transport format（XPORT）文件，需要把普通的SAS数据文件（.ssd和.sas7bdat）转换成Transport format(XPORT)文件，再用命令：read.xport()SPSS数据库：read.spss()可读入SPSS数据文件Epi info数据库：要给数据集一个名字，则是read.epiinfo("文件名.rec")->名称Stata数据库：R可读入Stata5,6,7的数据库读入数据文件后，使用数据集名$变量名，即可使用各个变量 >read.dta(“文件名.dta”) 读入数据文件后，使用数据集名$变量名，即可使用各个变量。>mean(data$age) 便是计算数据集 data中的变量age的均数。2.5 函数、循环与条件表达式2.5.1 编写函数句法是：函数名 = function (参数1，参数2…) { 函数体 函数返回值} 对于这类只有一个算术式的简单函数，也要不要{}>mean(data$age) 便是计算数据集 data中的变量age的均数。 若不使用圆括号，直接输入函数名，按回车键将显示函数的定义式：单参数：使函数个性化，可使用单参数，函数将会根据参数的不同，返回值不同> welcome.sb = function(names) print(paste("welcome",names,"to use R")) >welcome.sb("Mr fang") [1] "welcome Mr fang to use R" >welcome.sb("Mr Wang") [1] "welcome Mr Wang to use R" 默认参数：即不输入任何参数函数的默认参数> welcome.sb=function(names="Mr fang")print(paste("welcome", names,"to use R")) >welcome.sb() [1] "welcome Mr fang to use R" 当函数体的表达式超过一个时，要用｛｝封起来2.5.2 for循环for循环的句法是：for (变量 in取值向量) { 表达式… }
