感谢brain 老师和 所有跟制作教程有关 评论区热心解答问题 的人
results：输入的代码+结果 如果运行出错，也会显示
Command：输入代码 regrets y x 课程之后会用do file 可以替代 command
Review：回顾之前的操作 第一栏command 之前用过的代码
第二栏RC error code 如果之前的代码错了，RC 会是红色
Variables：变量
Properties：变量的类型是数值型、逻辑型
Note：还有注意的
Log：将所有命令和结果记录
Viewer：包含了所有的help file，
Graph：对绘制的图像进行更改
Do-file：可以将代码保存至此，具有可重复性
第三讲 初步数据观测
拿到数据集之后：
有多少个变量-列，观测值-行
符合条件的观测值有多少个？
一人一条数据？一人多条数据？
Describe：描述 count：计数 isid：是不是id unique：特殊的
选择内置数据集 1978 automobile data
[command] sysuse auto, clear
1.describe 种类（数值型），格式（保留几位小数），标签
syntax（语法）：describe [varlist] 变量名单=一个/多个变量
attention：[] 意味着 varlist 不是必须的，可以只输入 describe
如果只想看数据集的整体信息，不想看详细的 [command] describe, short
看每个变量的信息 [command] describe price
2.count 直接告诉该数据集有多少行rows observations
syntax：count
可以与logic 连用
syntax：count if
3.isid (is ID or not?) ID 是独特的，可以区分每一行观测值的；是否是这个数据的可识别标签
syntax：isid varlist 没有[ ] 因此 varlist是必须输入的内容
如果数据集，每人只有一个观测值，通过ID（每人的）就能区分每一行是哪个人的信息，若一人有两条或多条观测值，只看ID是区分不了一人的每一行数据，只能知道这些一行属于一个人。此时ID不是独特的，可以区分每一行数据的值
意思是有些观测值的mpg是相同的
mpg 和 weight相结合 是unique 的 identifier
也没有报错，证明大家的 price都是独特的
4.unique
盗版的 java installation not found 外部命令手动安装
网址“Index of /RePEc/bocode/u (bc.edu)”
ado pkg sthlp hlp 四种格式的后缀 如果有hlp 尽量把4个全部下载；没有hlp 前三个即可，【将链接另存为】
然后将【保存的链接】挪到【D:\stata17\ado\base\u】对应的首字母目录下，重启stata 即可
syntax：unique varlist
74 record 中 只有 21 value，所以 not uniquely identify，isid mpg就会报错
第四讲 统计描述指标1
上节课：变量数量、观测值数量、观测值筛选、具有独特值的数据量 describe、count、isid、unique
sysuse auto, clear
1.codebook 数据字典
一般拿到一个新的data set后，看感兴趣的变量的数据字典
overview of variable type, stats, number of missing/ unique values 类型，统计量，缺失值，特异值
syntax: codebook [varlist] [if] [in] [,options] [ ]内的不是必须输入项目
codebook该数据集中所有变量的信息
varlist：变量名单=一个或多个变量 if：逻辑判断 in：第几到第几个观测值 options：跟在逗号后面，一些可以自定义的选项
codebook price
codebook price if price >5000
codebook price in 10/20
in 10 看第10个观测值的codebook
in 10/30 …第10到30个…
in 10/l …第10到last（最后一个）观测值…
in f/10 …第1（first）到10…
help codebook
codebook, mv
不过一般codebook很少用到options，因为一般只是初步看看观测值
2.summarize
print summary statistics（mean，stdev，min，max）for variables
syntax：summarize [varlist] [if] [in] [weight] [, options]
推荐使用 sum 或 sum，完全等同于summarize
summarize ( sum / sum) price
detail 看一些额外的统计量，例子见下
sum price, detail （【detail】跟在【，】后面因为它是option）
|
偏度 丰度 |
|
最便宜的，第二便宜的… |
第五讲 统计描述指标2
histogram 直方图 graph box（箱形图；箱形图box-and-whisker plot） violin plot 小提琴图
sysuse auto, clear
1.histogram 直方图
varname：不是varlist，因此只能有一个变量
histogram可以简写成hist（观察下划线
syntax: histogram varname [if] [in] [weight] [, [continuous_opts | discrete_opts] options]
hist price
分成了8组，纵轴默认是density数据的密度，
hist price, freq hist price, percent hist price, frac
frequency 频数 percent 百分比 fraction 小数（也就是percent/100 总数是1，而不是100%）
可以改变bin的宽度（Bins 是您想要将所有数据划分成的间隔数，以便它可以在直方图上显示为条形。）
continuous 连续性变量，指定多少个bin，指定width宽度
discrete 离散型变量
hist price, freq bin(5) 设定5个bin，stata自动计算每个bin的宽度
还可以加一些密度曲线，如果样本量足够大，组距就会变得很短，阶梯折线变为曲线
hist price, freq bin(5) normal 拟合正态曲线
hist price, by(foreign) 按照foreign分组绘制
2.graph box 或者 graph hbox 箱线图（horizontal box横着的箱子）
syntax：graph box yvars [if] [in] [weight] [, options] graph hbox yvars [if] [in] [weight] [, options]
【复习】箱线图中，median 中位数，upper hinge是75百分位数，lower hinge是25百分位数
upper adjacent value为upper hinge+1.5*IQR（IQR=75百分位数-25百分位数）
lower adjacent value 为lower hinge-1.5*IQR
若数据落在upper/lower adjacent value之外，称之为outside values离群值
graph box price graph hbox price
graph box price, over(foreign) 通过foreign分层
3.vioplot 小提琴图（不是自带的 需要安装
vioplot price
|
75百分位有这么多观测值 |
|
25百分位有这么多观测值 |
vioplot price, over (foreign)
第六讲 散点图 scatter plots
调整stata的主题
graph query, schemes
推荐的主题 s1mono 干净整洁适合杂志发表
set scheme s1mono （stata重启前设置为s1mono主题
set scheme s1mono, perm （永久设置为s1mono主题 permanent 永久的
1.scatter plot散点图
[twoway] scatter varlist [if] [in] [weight] [, options]
twoway -- xy两个变量绘制图
varlist -- 最基础的形式 twoway scatter y x
进阶形式 twoway scatter y1 y2 y3 … x（x与不同y变量之间的关系
twoway scatter mpg weight
twoway scatter weight length price
改变散点图的形状、颜色、大小等等
twoway scatter mpg weight, msymbol( ) mcolor( ) msize( )
msymbol 改变形状（help symbolstyle
mcolor 改变颜色（help colorstyle
msize 改变大小（help markersizestyle
twoway scatter mpg weight, msymbol(D) mcolor(blue) msize(medium)
option：by 按照某分类变量分别绘制散点图
twoway scatter mpg weight, by(foreign)
想把domestic 和 foreign 整合在一张图上
twoway(scatter mpg weight if foreign==0)(scatter mpg weight if foreign==1),legend(label(1"Domestic")label(2"Foreign")) 报错 原因1、2后面没写空格
twoway (scatter mpg weight if foreign==0)(scatter mpg weight if foreign==1), legend(label(1 "Domestic") label(2 "Foreign"))
if foreign==0 -- 国内生产的车domestic
每个scatter 用括号括了起来，是一个单独的plot，都跟在twoway后面，合并到一起
legend 设置两个label，第一个scatter plot设置为domestic，第二个是foreign
第七讲 双变量图像2
sysuse uslifeexp, clear
help twoway 查看 stata可以绘制的双变量图像
twoway 命令入门
twoway plot 变量 [if] [in] [, twoway_options]
plot：选择图像的种类，比如 scatter, line, connected, area, bar等
变量：这里可以写一个或多个y变量，一个x变量。最后一个是x变量，之前的为y变量。
if：定义所取某一个自变量的范围
in：定义所取观测值的范围
twoway_options：可以定义图像的“美观”部分，例如坐标轴范围、标题、注释、标签等等
【上节课】
twoway scatter le year
twoway scatter le_male le_female year
le_male：y1 le_female：y2 year：x
折线图 twoway line y1 y2 … x
twoway line le year twoway line le_male le_female year
带数据标记的折线图 twoway connected y1 y2 … x
twoway connected le year twoway connected le_male le_female year
垂直线图 twoway dropline y1 y2 … x
twoway dropline le year twoway dropline le_male le_female year
脉冲图 twoway spike y1 y2 … x
twoway spike le year twoway spike le_male le_female year 图二
图二 twoway spike le_male le_female year
为什么只显示了female（灰色）？ twoway plot y1 y2 x，y2会覆盖在y1之上
twoway spike le_female le_male year 处理方式，将顺序调换
twoway spike le_female le_male year, scheme(s2color)
面积图 twoway area y1 y2 … x
twoway area le year twoway area le_female le_male year
LOWESS图（LOESS图）twoway lowess y x
twoway lowess le year lowess le year
第八讲 统计描述指标3
1.ci mean：连续变量mean（平均值）的置信区间
compute standard errors and confidence intervals
ci means [varlist] [if] [in] [weight] [, options]
cii means #obs#mean#sd [, level(#)]
默认95%置信区间（可通过 level(99)等进行更改）
ci mean mpg price, level(95)
cii mean 166 19509 4379, level(95)
直接告诉stata 观测值的数量=166 # 平均值=19509 # 标准差=4379standard deviation
2.proportion（可简写为prop）：分类变量mean的置信区间
分类变量的置信区间
方法1 Choice 1：ci proportions [varlist] [if] [in] [weight] [, prop_options options]
ci prop foreign
缺点：ci prop只能用于binary二分类变量
codebook rep78
value 1 2 3 4 5 ·(缺失值)，有这么多种分类，所以不是二分类
ci prop rep78
方法2 Choice：proportion分类变量的置信区间
proportion varlist [if] [in] [weight] [, options]
默认95%置信区间
prop foreign
二分类变量，所以proportion 加起来=1
prop rep78
number of observations = 69，去除了缺失值后的数量
proportion1 2 3 4 5加起来=1
prop foreign rep78
prop foreign prop foreign rep78
观测值数目=74 观测值数目=69
Domestic proportion =.6956522 Domestic proportion =.7027027
为什么会发生变化？
74-69=5，有5个缺失值
执行 prop foreign rep78【prop多个变量时，stata默认去除缺失值之后再计算proportion
【解决办法】prop foreign rep78, miss
但是！！报错了 没有出现讲课截图中的样子（下图）
3.pwcorr：变量的配对相关性
pwcorr [varlist] [if] [in] [weight] [, pwcorr_options]
pwcorr price headroom mpg displacement
对角线上的数据都是1，因为price-price
price-headroom：0.1145 正相关 price-mpg：-0.4686负相关
pwcorr price headroom mpg displacement, sig【sig 显著性 即p值，在每个相关性的下方显示出来】
不想看所有的相关性的p值，只想看p值<0.05的。
pwcorr price headroom mpg displacement, star(0.05)
加*的为p值<0.05者
4.graph matrix：绘制相关性矩阵（图像）
graph matrix varlist [if] [in] [weight] [, options]
graph matrix price headroom mpg displacement
|
可以观察到 对角线两侧是完全对称的，因此我们可以只看一半
graph matrix price headroom mpg displacement, half
第九讲 可重复性：引入 reproducibility
Log File 所有显示在results上的结果都保存为文档形式
类似于R markdown文件
Do file 合作者、提醒自己、存档
Do flie写一个read me file 需要安装fre
第十讲 可重复性：Log文件
Log command
log using filename [, append replace [ text | smcl ] name( logname) ]
filename：Log file 的名字
任选其一 append：如果文件存在，附加append在文件上
replace：如果文件存在，替换replace该文件
如果文件不存在，append和replace都会创造新的文件
log using yikahui, append
smcl：stata默认的log file格式，可以保存各种颜色
text：单色，便于在文本编辑器中打开
作者偏向于text，比较干净，体积较小
name(logname) 通常在一个stata运行的过程中，只有一个log file
已经有一个log file在运行了
有时，想创建不同log file给不同的合作者，又想在一个代码中完成并记录
因此：如果给log file起一个名字name( logname)，那么可以同时打开好几个log file
不同的log file给不同的合作者
log using Brian, name( Log_For_Brian)
log using Allen, name( Log_For_Allen)
log using Jimmy, name( Log_For_Jimmy)
于是打开了3个不同的log file
关闭log file
没有起名字（最常见）
log close
如果起了名字，需要指定名字
log close logname
log close Log_For_Brian
一口气关闭所有log file，包括起了名字的
log close _all
第十一讲 可重复性：Do file
打开文件 Lecture11.do
绿色文字为注释，在stata中并不会运行，可以协商【该command是做什么用的、提醒自己注意的事项】
*只能作为一行的标记，/* */ 可以跨越许多行，建议在每个do file的开头做一个【Purpose、 Author、 Date created、 Last modified】，告诉自己该 Do file的用处
选中 clear all（清除掉stata内存中的所有数据），点击do
打开 data browse，即发现是空的，代表stata内存里没有数据
建议创建一个文件夹，里面保存所有数据log file和do file，便于数据的可重复性
pwd 查看现在的工作路径
如果想改变工作路径 cd “D:\stata17\外部命令“
定位好工作路径之后，stata便会在该路径下寻找操作所需的data，
log using "Lecture11", replace
use auto.dta, clear 调用一个 .dta格式的数据集
（该图片可以看出 工作路径 的意义）
如果数据集是其他格式的，如excel、csc等，之后的课程会讲解如何进行数据导入
, clear 是一个比较保险的option，告诉stata：如果内存里已有数据，清除之后再导入
除了使用command，还可以使用窗口菜单，见下图
通过窗口菜单使用的command，可以放入 do file中：全选 复制 粘贴 do
do file中，一行command可以只选中一部分，就运行do
可以不选择任何一个command，直接运行该do file中全部command
第十二讲 t检验
t检验包括 单样本t检验 One-sample t test
独立样本t检验 Two-sample t test
配对样本t检验 Paired t test
因为是盗版的，不能用 webuse
想看price 的 mean 是否等于6000，如下图设置（by/if/in 设置各种条件判断）
单样本t检验 One-sample t test
ttest varname == # [if] [in] [, level (#)]
ttest price == 6000
level 默认95%水平
ttest price == 6000, level(99)
ttest price == 6000 if foreign == 0, level(90)
独立样本t检验 Two-sample t test
EG1：
EG2：
webuse fuel.dta, clear（因为是盗版 没法使用联网数据库
在收藏夹里有stata官网，检索相应的名字，即可找到.dta
>0.05不拒绝原假设
配对样本t检验 Paired t test
双样本t检验
一个变量分组比较 ttest varname [if] [in], by(groupvar) [level(#)]
例子：ttest price, by(foreign)
两个变量比较
非配对样本：ttest varname1 == varname2 [if] [in], unpaired [level(#)]
例子：ttest mpg1 == mpg2, unpaired
配对样本：ttest varname1 == vaarname2 [if] [in], [level(#)]
例子：ttest mpg1 == mpg2
第十三讲 卡方检验/Fisher精确检验
sysuse nlsw88, clear
卡方检验 Chi-squared test
Pr<0.005 不同婚姻的race分布不同；不同race的婚姻状况不同
frequency 487…
Chi-square contribution 16.7 9.3 …
expected frequency期望频数（有时又叫理论频数 586 1051…
在 key 可以看排列的内容种类
Fisher精确检验 Fisher’s exact test
通常，当列联表中理论频数（期望频数）<5时，我们可以增加样本量、删去理论频数太少的行或列、或者合并某些行或列；也可以使用Fisher精确检验
任何样本量都可以使用Fisher精确检验
使用代码进行卡方检验、Fisher精确检验
卡方检验：tabulate var1 var2, chi2
每个单元格对卡方检验的贡献：tabulate var1 var2, cchi2 chi2 （cell chi-square==cchi cchi 和chi 顺序可以互换）
理论频数（期望频数）：tabulate var1 var2, chi2 expected
Fisher精确检验：tabulate var1 var2, exact
注意：
1 总例数≥40，所有理论频数≥5，看Pearson Chi-Square结果
2 总例数≥40，出现1个理论频数≥1且<5，X2检验需进行连续性校正，以Continuity Correction结果为准
3 总例数≥40，至少2个理论频数≥1且<5，看Fisher’s Exact Test结果
4 总例数<40或者出现理论频数<1，看Fisher’s Exact Test结果
为什么可以不用连续性校正？
- stata不自带卡方检验的连续性校正
- stata有用户自写的package，但不推荐；卡方检验的连续性校正并不是必须的，也不是最推荐的方法
- 在样本量足够大的时候，使用卡方检验时，是否使用卡方检验的连续性校正区别很小；使用Fisher精确检验也是没有问题的
- 在样本量小的时候（通常是某个格子期望频数<5），可直接使用Fisher精确检验，亦不需要使用“卡方检验+连续性校正”
第十四讲 RR值计算
同样 去stata database 中检索 csxmpl.dta 下载（如遇安全问题 选择 “保留”）
sysuse csxmpl, clear
list
数据的结构
具有三个变量：case：1==case 0==non-case
exp：1==exposed 0==non-exposed 暴露/非暴露
pop 人数population
选择pop人数进行加权，否则 stata四行每一行都会当作观测值，但我们知道每一行是许多观测值的合并
fweight 频数加权
risk difference（相减）与0参照 置信区间 -0.7~-0.1 不包括0 所以p<0.05
risk ratio（相除）与1参照 置信区间0.2~0.9 不包括1 所以p<0.05
队列研究计算器
使用代码进行RR值的计算
cs var_case var_exposed [if] [in] [weight] [, cs_options]
cs case exp [fweight = pop]
csi #a #b #c #d [, csi_options] 
zhuyi注意 #不输入
第十五讲 OR值计算
OR值和RR值的区别？
RR值：Cohort study或者RCT中，研究者前瞻性地观察“暴露组”和“非暴露组”的发病情况，之后通过RR来评价“暴露组”研究对象的发病风险是“非暴露组”研究对象的多少倍？
OR值：在回顾性研究（如case-control）中，研究对象是已经患病的“病例组”和未患病的“对照组”，研究者回顾性地调查病例组和对照组的暴露情况，因此无法计算发病率等指标。
最终目的：知道相对风险
OR是否可以估计RR值？
当终点事件发生率较低时，OR可以近似为RR（<15% 不是金标准）
当终点事件发生率较高时，OR会“夸大”RR“值
OR值相对于RR值“更远离1“
当RR值大于1时，OR大于RR（1<RR<OR）
当RR值小于1时，OR小于RR（OR<RR<1）
终点事件发生率越高时，OR越会overestimate
对于队列研究/RCT，可以报告OR么？
可以，但是：RR对于效应值的估计更加准确，对于临床意义的解释更加明确（病人理解，吃药 比 不吃药 使疾病风险降低50% √ 吃药比值比是50% ×）
regression model回归模型 中：对于结局是二分类变量的研究，logistic回归只能提供OR值，不能提供RR值（之后会讲：当结局发生率高时，应该使用log-binomial回归或者使用带有稳健方差估计的泊松回归，直接提供RR值
接上节课RR内容
如果看risk ratio 该暴露降低了49%的风险，odds ratio 降低了87%，因此odds ratio是夸大了暴露因素的影响
强调：能提供RR值时（队列研究、RCT）请不要提供OR值
回顾结束，本节课内容————————————————————
调入并观察数据（汇总数据）
webuse ccxmpl, clear （依旧是手动 sysuse ccxmpl, clear）
对话框操作（使用数据集算OR）
病例对照优势比计算器
使用代码进行OR值的计算
cc var_case var_exposed [if] [in] [weight] [, cc_options]
cc case exp [fweight = pop]
cci #a #b #c #d [, cci_options] 
第十六讲 方差分析ANOVA
单因素方差分析 F值 与 单因素回归分析 t值 是相同的 t2即F，p值也一样
多因素方差分析 不如直接使用 多因素回归分析 后者能知道每一个变量对结局的贡献
回归分析中 能更加方便的看两个变量是否存在交互作用 ANCOVA（co-variance
下节课：线性回归 logistic回归 log-binomial回归
方差分析的假设
假设1 y变量为连续变量
2 有一个包含2个及以上分类、且组别间相互独立的x变量
3 每组间和组内的观测值相互独立（123在实验设计阶段尽力）
4 每组内没有明显异常值
5 每组内y变量符合正态分布
6 进行方差齐性检验，观察每组的方差是否相等
sysuse systolic, clear
假设4：每组内有无明显异常值？
方法：绘制 Boxplot（violin plot也可以
假设5：每组内y变量符号正态分布？
点击 data browser看一下数据结构
systolic作为结局变量 drug分1234
codebook drug
sum systolic, detail
假设4 验证
increment in systolic可能为负 第三个的下限能包括进该值，故暂不认为是异常值
假设5
- 偏度和峰度正态性检验
p值 0.1331 > 0.05 不能拒绝 符合正态分布 这一原假设
②Swilk-Shapiro-Wilk正态性检验
p值 0.37331 > 0.05 不能拒绝 符合正态分布 这一原假设
③sfrancia- Shapiro-Francia正态性检验
同上
三种检测方式任选其一即可 很多人使用②，提前说好用哪一个，
完成假设4 5 的检验之后，可以开始oneway ANOVA
df 3 分子的自由度 54 分母的自由度 F 3 54分别是多少？肯定远远小于9.09 因为p 0.0001
所以四个组里至少有两个组的平均值是有显著差异的
假设6 方差齐性检验 观察每组方差是否相等
Bartlett’s test for equal variances：p值0.800 >0.05 通过方差齐性检验，因此四个组方差齐
到底是哪两个组呢？使用两两比较 示范用的Bonferroni
p值1.000 已经经过Bonferroni 处理 所以p值跟0.05比较
1-3 2-3 4-1 4-2 的p值<0.05 都是有显著差别
row mean- col mean 组2-组1=-053 组3-组2=-17.31 等等
使用代码
oneway systolic drug oneway y x 与之后的 regression 相似 regress y x regress y x1 x2
oneway systolic drug, bonferroni
oneway systolic drug, bonferroni tabulate
twoway方差分析
第十七讲 简单线性回归 simple linear regression
线性回归的假设
假设1：y是连续变量
2：x可以被定义为连续变量（这一条不严格 x可以是哑变量
3：y和x之间存在线性关系
4：具有相互独立的观测值
5：不存在显著的outlier
6：等方差性
7：residual近似正态分布
假设3：因变量和自变量之间存在线性关系
散点图/Lowess图（分别对应 第六讲：散点图 第七讲：双变量图像2）
如果不符合线性关系怎么办？
Spline：一次方向，分段fit直线
Quadratic：二次方项
Cubic：三次方项
Restricted cubic：三次方项（头、尾近乎直线）
…
假设4：具有相互独立的观测值
可以使用杜宾-瓦特森（Durbin-Watson）统计量
stata对于非time series（时间序列）数据不设有这个统计量的检测
更多情况下，不需要测量这个假设是否成立（实验设计阶段 即完成
如果是相互关联的观测值：GEE模型，Multi-level模型
假设5：不存在显著的异常值
Boxplot或者violin plot（第五讲
假设6：等方差性
使用Residual-versus-fitted plot（纵坐标residual横坐标fitted
代码：rvfplot
窗口菜单：Statistics → Postestimation
假设7：回归残差近似正态分布
step1：得到残差
step2：使用正态分布的检验方法
直方图
使用qq plot观测
偏度峰度、Shapiro-Wilk检验、Shapiro-Francia检验
data browser
假设3：因变量和自变量之间存在线性关系
set scheme s1mono
twoway scatter mpg weight
基本是一个线性关系
twoway lowess mpg weight
lowess mpg weight 将lowess 和 散点图 显示在一起
假设4：具有相互独立的观测值（假定符合，因为都是不同的车嘛
假设5：不存在显著的异常值
Boxplot或者violin plot（第五讲
最上面的横线：75分位+IQR 并没有多出去多少，所以姑且不算异常值
假设6：等方差性
regress mpg weight 建议回归分析直接敲代码 因为简单：regress y x / regress y x1 x2…
F（1，72）整个模型有无共线：比如有好几个x，测所有变量的β值都等于0，（？？）
R-squared纳入模型的变量能够解释y变量的多少 65% 64%
当纳入的变量越多，R-squared↑，但大多时候这对模型没有意义，因为可能有共线性、过拟合等等；Adjusted R-squared 进行过校正更好显示贡献
截距项39.44（weight=0时 mpg=39.44）
rvfplot
fitted value 即 predicted y value
看residual是否分布在0两侧
假设7：回归残差近似正态分布
predict cancha, residual
cancha 随便取一个新变量的名字（cancha残差
data browser
做一个cancha残差的直方图
jinsi
近似正态分布
也可以用qq plot观测
如果贴合这条直线就是完全正态分布，
使用偏度峰度，Shapiro-Wilk检验、Shapiro-Francia检验
p值<0.05
第十八讲 多元线性回归
sysuse auto.dta, clear
price=β0+β1weight regress price weight
β1：汽车的重量每增加一个单位，售价增加2.04个单位（95%CI：1.29，2.80）
p=0.000？× p<0.001 √
β0怎么解释？：当汽车的重量=0时 汽车的售价=-6.70，但这不合理：汽车重量不可能=0，售价也不可能为负数
price=β0+β1weight+β2length regress price weight length
β1：在控制了汽车的长度后，汽车的重量每增加一个单位，售价增加4.70个单位（95%CI：2.46，6.94）
β2：…重量…，长度…，售价降低97.96个单位（95%CI：-176.08，-19.85）
price=β0+β1weight+β2length+β3mpg regress price weight length mpg
汽车的价格是否和修理次数有关？
codebook rep78
5个缺失值
price=β0+β1weight+β2length+β3mpg+β4rep78 regress price weight length mpg rep78
β4：在控制了汽车的重量、长度、里程后，汽车的修理次数每增加一个单位，售价增加910.99个单位（95%CI：302.62，1519.35）
这里rep78当作了连续变量，但应该当作多分类变量--
Dummy variable：哑变量
当自变量x为多分类变量时，例如职业、学历、血型、疾病严重程度等等，此时仅用一个回归系数来解释多分类变量之间的变化关系，及其对因变量y的影响，就显得不太理想
在多分类变量前加上“i.“，告诉Stata这不是连续变量
regress price weight length mpg i.rep78
rep78=1 被当作对照组，
Irep78_2：在控制了汽车的重量、长度、里程后，汽车的修理次数两次和一次相比，售价增加1137.28个单位（95%CI：-2468.70，4747.27）
回归分析中需要注意的事情！
1、时刻注意纳入模型的观测值数量
rep78中有5个缺失值，
Case-wise Deletion 只用纳入模型的变量都没有缺失值，这个观测值才会纳入回归模型分析！
也就是 该observation 在 price/ weight/ length/ mpg/ rep78五个变量中都没有缺失值，才能纳入
2、如何让β0有意义？
汽车平均重量：3019.459（sum weight）
转化为 price=β0+β1（weight-3019.459）
β0：当weight=平均重量时，汽车的售价
gen weight_center=weight-3019.459
regress price weight_center
窗口菜单操作
第十九讲 二分类Logistic回归 Logistic regression(binary outcome)
回顾第十五讲 OR是否可以估计RR值？“当终点事件发生率较高时，OR会夸大RR值”
横坐标 终点事件发生率 纵坐标 OR值
可以看到，终点事件发生率20+时，RR值为3（真实的relative risk），OR值为8，夸大了
（偏离RR值越大 不论是减少 或增加 都称之为overestimate，因为相关性 也有正负 值的绝对值越大 相关性越强）
当outcome发生率>15%时，logistic regression得出的OR值会overestimate实际的RR值
解决办法
- 传统的Logistic regression（得出OR值）
- Mantel-Haenszel方法（得出RR值）
- Poisson regression with robust variance estimate（泊松回归
- 新方法：
做Logistic回归之前检查
条件1：Y是二分类分类变量
条件2：Y的发生率<15%
（条件1、2就是区分 Logistic 和 log binominal 的区别
使用二分类Logistic模型前，需判断是否满足以下7项假设
假设1：因变量（结局) 是二分类变量
2：有至少1个自变量，自变量可以是连续变量，也可以是分类变量。
3：每条观测间相互独立。分类变量（包括因变量和自变量) 的分类必须全面 且 每一个分类间互斥。
4: 最小样本量要求为自变量数目的15倍，但一些研究者认为样本量应达到自变量数目的50倍。
5：连续的自变量与因变量的logit转换值之间存在线性关系。
6：自变量之间无多重共线性。
7：没有明显的离群点、杠杆点和强影响点。
sysuse lbw, clear
变量low（<2500g）是结局事件，想了解什么因素和孩子的low birthweight 相关
codebook low
tab low
发生率31.22%>15%
disclaimer免责声明：本次课的数据集中，结局事件的发生率远大于15%，应使用Log-binomial模型进行分析。本次课使用Logistic regression进行分析仅为了讲解如何使用stata进行操作，以及为下节课讲解Log-binomial进行铺垫
Logistic regression
Stata code：logistic y x1 x2 x3…
窗口菜单：
Model1：low=β0+β1age（母亲age与孩子low weight是否有关
logistic low age
β1：母亲年龄每增加1岁，孩子低体重风险是之前的0.95倍（95%CI：0.89，1.01）
Model2：low=β0+β1age+β2smoke
logistic low age i.smoke
β1：控制了母亲吸烟状况后，母亲年龄每增加1岁，孩子低体重风险是之前的0.95倍（95%CI：0.89，1.01）
Model3：low=β0+β1age+β2smoke+β3race
logistic low age i.smoke i.race
β1：控制了母亲吸烟状况和种族后，母亲年龄每增加1岁，孩子低体重风险是之前的0.97倍（95%CI：0.90，1.03）
β2：控制age和race后，吸烟者与不吸烟者相比，孩子低体重风险是3.01倍（95%CI：1.45，6.23）
β3：控制age和smoke后，black与white相比，孩子低体重风险是2.74倍（95%CI：1.04，7.23）
【stata会默认将值最小的作为多分类变量中的reference，此处即race中的white】
窗口操作
第二十讲 广义线性模型 generalized linear model, GLM
回顾：大于15%会产生bias，这个问题只出现在我们想用odds去estimate risk的时候才会发生
- 在cohort study和clinical trail中，outcome是结局事件的发生率
- 在cross-sectional中，outcome是结局事件的prevalence
- 在case-control中，我们不报告risk，只报告odds，因此就没有这个问题
在一个研究中，我们只想报告odds，而不去estimate risk，无论结局事件发生率是多少，logistic回归都是valid。
线性回归（lecture17-18）：Expected Y=β0+β1*X1+…+βn*Xn（不用对Y进行转换
广义线性模型：f（expected Y）=β0+β1*X1+…+βn*Xn（f 意味着function）
- linear regression is a special case of GLM
- same for logistic regression, Poisson regression, Log-Binomial regression
- even for Cox regression, GEE model
- 在GLM中，我们转化Y，使得转化后的Y和X们呈线性关系
simple GLM 与 Multiple GLM
简单广义线性模型 f(expected Y)=β0+β1*X1
多元广义线性模型 f(expected Y)= β0+β1*X1+…+βn*Xn
x变量可以是各种类型的变量（连续、分类）
in GLM，we have to specify 2 things
family == Y变量如何分布？ 比如正态分布、泊松分布、二项分布
link ==把Y怎么转化才能和X成linear关系？
Linear regression 最简单的线性回归
Y variable ①continuous variable
②如何分布？ 高斯分布（Gaussian Distribution）
它的family是：我们有一个assumption, 在线性回归中Y变量
应该是正态分布或高斯分布
model ①Expected Y=β0+β1*X1+…+βn*Xn
②Y怎么转化才和X成linear关系？不用转化！link: Identity(恒等 意思：不转化它便相等，即为恒等（？)
Logistic regression
Y variable ①Binary variable
②如何分布？ 二项分布（Binomial）
model ①ln [ P ( Y=1 )/ P ( Y=0 ) ]=β0+β1*X1+…+βn*Xn （ P(Y=0)可以写作 1-P( Y=1 ) 因为是二项分布）
②Y怎么转化才和X成linear关系？link: Logit
Poisson regression
Y variable ①count variable：整数（1, 2, 3…n）
②如何分布？泊松分布（Poisson）
model ①ln [risk of event]= β0+β1*X1+…+βn*Xn
②Y怎么转化才和X成linear关系？link: Log
Log-binomial regression
Y variable ①binary variable
②如何分布？二项分布（Binomial）
model ①ln [ P ( Y=1) ]= β0+β1*X1+…+βn*Xn
②Y怎么转化才和X成linear关系？link: Log
help glm 可以查看该表内容
Linear regression: lecture 17-18 regress y x1 x2 x3
Logistic regression: lecture 19 logistic y x1 x2 x3
Poisson regression 右侧二维码 poisson y x1 x2 x3
glm command in stata 就像一把雨伞 什么都可以装在里面
umbrella command
linear regression: family ( gaussian ) link( identity)
logistic regression: family ( binomial ) link( logit)
poisson regression: family ( poisson ) link( log )
log-binomial regression: family ( binomial) link( log )
Linear regression
regress price weight length mpg i.rep78
glm price weight length mpg i.rep78, family(gaussian) link(identity)
Logistic regression
logistic low age i.smoke i.race
glm low age i.smoke i.race, family( binomial) link( logit) 【logistic regression 报告odds ratio，e^n 报告的是n】
glm low age i.smoke i.race, family( binomial) link( logit) eform 【但我们想知道odds ratio 即e^n ，本行command报告的是e^n】
对比
regress price weight length mpg i.rep78 (见lecture 18)
glm price weight length mpg i.rep78, family(gaussian) link(identity)
logistic low age i.smoke i.race (见 lecture 19)
glm low age i.smoke i.race, family( binomial) link( logit) eform
Log-binomial Model
glm low age i.smoke i.race, family( binomial) link( log)
alternative to Logistic regression【是logistic regression 的替代
问题：fail to converge (对于一个command 一直跑代码 停不下来
解决办法solution：Poisson regression with robust variance estimate 使用带有稳健方差估计的泊松回归
glm low age i.smoke i.race, family( poisson) link( log) robust ×
glm low age i.smoke i.race,family(poisson) link(log) robust eform √
【此处老师讲的代码和图不匹配 应该是eform这一行 对应视频中的IRR的图。感谢评论区
【incidence rate ratio 即IRR，此处当成relative risk也可以
窗口菜单
第二十一讲 Kaplan-Meier曲线和Log-Rank检验
回顾：Survival analysis in use 做生存分析时的目标：
①描述一个组内个体的生存时间：寿命表法（life tables methods）、非参数Kaplan-Meier曲线（lecture 21）（后者逐渐成为主流）
②比较两个或多个组的生存时间：Log-rank test（lecture 21）
③研究生存时间和变量之间的关系：半参数Cox比例风险模型（lecture 22）、参数生存分析模型（too advanced, not covered in these 30 lectures）
K-M曲线的历史：
介绍了一种全新的、解决随访期间Right Censoring的问题的生存分析方法
特点：精确地记录并利用每个个体发生终点事件的具体时间，在任何一个终点事件发生的时间点计算出一个新的、基于之前所有信息的Cumulative survival
优点：①相比于life-table method，更加充分地利用了信息，给出更准确的统计量
②非参数估计方法：不要求总体的分布形式，因此非常适合（初期的）生存分析时使用
③K-M曲线可以很直观的表现出两组或多组的生存率或死亡率，适合在文章中展示
手动下载drugtr.dta
导入数据：sysuse drugtr, clear
恢复成为普通的数据格式：stset, clear
stata已经将这个数据集设置成了生存数据的格式，导入数据集后，将数据集恢复成普通的数据格式，才是我们在临床研究中见到的数据结构
step1 数据集的初步观测
list values of variables
list
list in 5/10 呈现5 6 7 8 9 10这六个观测值
有四个变量，每个的观测值可以看见
codebook drug
0=placebo 安慰剂 所以1=试验药
codebook studytime
step2 声明数据为生存分析数据代码操作
告诉stata：终点事件（failure variable），随访时间（time variable）
stset timevar, failure ( failvar [==numlist ] )
timevar：随访时间变量
failvar：终点事件变量
numlist：终点事件变量中，哪些值代表发生了终点事件（例如1代表终点事件，numlist即为1）
stset studytime, failure ( died==1)
窗口菜单
多记录的ID变量multiple-record ID variance：
- 若不指定ID，stata则认为每条观测值代表一位患者（临床研究常见如此
- 若一位患者有time-varying exposure，则其有多条观测值，应在这栏中告诉stata用哪个变量区分不同患者（通常是代表患者ID的变量）
- 这种情况更常见于队列研究中，而在实际的临床分析中并不常见（因为患者不会一会儿吃placebo 一会儿吃实验药）
step3 生存数据再观测
注：必须要在指定数据集为生存分析数据集之后（stset之后）才能使用任何其他的st开始的命令。
stsum
50%中位生存时间，很在乎
窗口菜单
stdescribe
窗口操作
sts list
stata会在每一个时间点上计算出一个新的tumor list survival
窗口操作
step4 K-M曲线的绘制
可以画生存曲线 或者 死亡曲线
代码：sts graph, by( var )
sts graph, by( drug )
sts graph if age<50, by(drug) 绘制age<50的人
sts graph, by(drug) risktable
期待：高级绘图1-2 lecture29、30 讲解
fancy 的图捏！
step4 检验组间差别（Log-Rank Test）
sts test drug 比较drug=1 和 drug=0这两个组的survivor functions是否一致
p<0.0001 因此有显著差别，得出结论：实验药能够显著提高生存率
第二十二讲 Cox回归和比例风险假定检验
回顾：Kaplan-Meier曲线： 描述一个组内个体的生存时间
对于两组或多组患者的生存率、死亡率进行直观的比较
Log-rank检验：比较两个或多个组的生存时间
得到统计学上的相应证据，即p值
但：Kaplan-Meier曲线是一种非参数的估计方法，不能控制潜在的混杂因素对于结局事件的影响，也不能看出某种因素，如患者的年龄、性别、种族等，是否是结局发生的独立因素或危险因素等，
question：在一个抗癌药物的Clinical Trial中，48名患者被随机分配到新药组（28人）和安慰剂组（20人），研究人员想知道新药是否影响患者的生存情况
sysuse drugtr, clear
提示：数据初步观测 和 转化为生存数据格式 在上节课已讲，在此不赘述
step1 Cox回归
窗口菜单
statistics → survival analysis → regression models → Cox proportional hazards model
我们将数据转化为survival data时（stset），指定过终点事件（failure variable）、时间变量（time variable），我们在此只需告诉stata 放入哪些x变量
drug的coefficient 新药组终点事件发生风险是安慰剂组的13.3%（95%CI：5.6%，31.4%）
“与安慰剂组相比，新药组发生终点事件的风险降低了86.7%”
drug的coefficient 在控制了患者年龄age后，新药组终点事件发生风险是安慰剂组的10.5%（95%CI：4.3%，25.6%）
age的coefficient 在控制了治疗方法drug后，患者年龄每增加1岁，发生终点事件的风险增加12.0%（95%CI：4.1%，20.5%）
代码操作
stcox va1 var2 var3 … [if] [in] [, options]
再提醒：①必须指定data为survival data后（stset）才能使用任何st开头的命令
②因一开始已将data转化为survival data时，指定过终点事件failure variable、时间变量time variable，我们在此只需设置回归方程中independent variable即可
example：stcox drug age
stcox drug age if age<50
step2 PH假定检验
窗口菜单
statics → survival analysis → regression models → test proportional-hazards assumption
PH是一个必须要紧跟在回归模型之后的，一个对于模型的评价和估计
其零假设：纳入Cox回归模型的变量满足PH假定
p=0.8064>0.05，不能拒绝零假设，所以满足PH
注：若前面没有进行Cox回归，输入代码 estat phtest，会报错
通过图像观测变量是否满足PH假定
statics → survival analysis → regression models → graphically assess proportional-hazards assumption
set scheme s1mono
如果两条线近乎平行，则满足PH假定
调整了age之后，变化不大，依旧近乎平行
代码操作
eatat phtest
再强调：该命令需紧跟在Cox回归的命令之后，否则stata不知道检验哪个回归的PH假定
stphplot, by ( var1 ) adjust ( var3 var3… )
var1是自变量名，var2是希望控制的变量。注：这个命令不一定要跟在cox回归之后
stphplot, by ( drug ) adjust ( age )
如果不满足PH假定怎么办？
- 一般只要两组生存曲线趋势一致、不明显交叉即可判定PH假定成立
- 如果PH假定不成立，可以加上时间（time）和暴露（exposure，
比如本例之中的drug）的交互项（interaction term），即time*exposure
【生成一个新变量放入模型，观察药物在不同时间的作用效果放松PH假定
- 我们也可以对于不同的时间段分别分析（e.g. 0-10, 10-20, >20）
- 参数生存分析模型：streg进行参数生存分析
第二十三讲 广义估计方程 generalized estimating equation GEE
examples of longitudinal data analysis LDA，纵向数据分析
child BP measured at each annual visit from 3 to 9 years old
infant weight measured at 3,6,9,12
in a 3-arm cross-over trial, short chain fatty acid measured at the end of each other
(cross-over trial 理解：实验中三个组 高糖高脂高蛋白，先把人进组，吃一段时间高脂，洗脱，吃一段时间高糖，洗脱吃一段时间高蛋白)
纵向数据分析到底比横断面分析强在哪儿？
a classical example: Association between age and reading ability? 孩子的年龄与阅读理解能力
若是横断面研究，看起来似乎年龄↑阅读能力↓
不论基线水平，对一个孩子而言，年龄↑阅读能力↑
更有说服力，虽然和横断面一致
advantage of LDA 优势
establish temporal trends (时序：A在B之前发生)
separate cohort effects from aging effects
children have different baseline values in reading abilities
the trajectories of reading abilities differs by people
from mathematical perspective
若不视作一个个体的mpg1 mpg2，而分开
p值变小了，差值diff没变，故是standard error变小了（1.23→0.78），
Var ( X-Y ) = Var ( 1 x X +(-1) Y )=( 12 ) Var ( X ) + ( -12 ) Var ( Y ) + 2 ( 1 ) ( -1 ) Cov (X, Y)
非配对t检验中，红色=0，配对t检验中，红色≠0，因此Var↑
LDA中，when estimating the “change”, taking into account the correlation between observations will give us small SE ( standard error )
线性回归的假设中，观测值需要有相互独立（lecture17 假设4）
within an individual, the observations are correlated
the next question is: how are they correlated?（观测值如何关联）
independent correlation structure 最常用
如果实际情况不符合independent correlation structure，Robust variance estimate if worried about mis-specification of the correlation structure
use QIC(an alternative for AIC for GEE models) to see which model works best
data: NLS women 14-24 in 1948 (national longitudinal surveys)
sysuse union, clear
观察数据list in 1/10
ID=1的人在72年随访，在71年没有，id=2在71年随访
告诉stata：I want to use GEE!
xtset studyid timevar
studyid: unique ID for each participant
timevar: timevar will usually be a variable that counts 1, 2, …, and is to be interpreted as first year of survey, second year, …, or first month of treatment, second month.
xtset id year
GEE command
xtgee y x1 x2 x3 … [if] [in] [weight], family (family) link (link) corr(correlation structure) robust
family: lecture 20 ; or help xtgee
link: lecture 20 ; or help xtgee
correlation structure
xtgee union age grade not_smsa south, family(binomial)link (logit) corr(ind) robust
第二十四讲 有序多分类Logistic回归 logistic regression(ordinal)
多分类变量分为有序多分类和无序多分类变量
有序多分类变量：疾病分期（1-4期）；疼痛情况（1-10：1最轻，10最痛）
无序多分类变量：居住地（1东/2南/3西/4北）；手机品牌（1苹果/2三星/3华为/4小米/5其他）
有序多分类的原理
将y变量的n个分类拆成n-1个二分类Logistic回归
今天例子中的y变量有5个level：Excellent; good; average; fair; poor拆分成：
poor(1) vs excellent + good + average + fair(0)
fair + poor(1) vs excellent + good + average(0)
average + fair + poor(1) vs excellent + good(0)
good + average + fair + poor(1) vs excellent(0)
proportional odds假定
- 多个二元logistic回归中，除了β0以外的系数相等 即OR1=OR2=OR3=OR4
odds ( poor) /odds ( excellent + good + average + fair ) OR1
odds ( fair + poor) /odds ( excellent + good + average ) OR2
odds ( average + fair + poor ) /odds ( excellent + good) OR3
odds ( good + average + fair + poor) /odds ( excellent ) OR4
- proportional odds假定是否成立更多是由研究问题的自身性质决定，可以用数据进行检测（下节课），但数据本身可能有Bias
- 如果该假定不成立：当作无序多分类logistic回归（下节课）
本节课使用数据 1977年汽车修理记录数据
sysuse fullauto, clear
codebook rep77 outcome：汽车维修状况
结局变量5个level
codebook foreign exposure：是否为进口车
step1 Chi-squared test(卡方检验)
H0：车辆是否为进口车和车辆维修状况没有关系
tab foreign rep77, chi2
p=0.008<0.05 拒绝零假设，得出结论：车辆是否为进口车和车辆维修状况有关系
→下一步：如何量化这种关系？
step2 Ordinal Logistic Regression
ologit y x1 x2 x3 [if] [in] [weight] [, options]
常用的options包括or（若不加or，stata给出结果为coefficient，得不到OR
examples ologit rep77 foreign
ologit rep77 foreign, or
ologit rep77 foreign length mpg, or
ologit rep77 foreign
1.46 进口车（foreign=1）和国产车（foreign=0）比：
odds (poor) / odds (excellent + good + average + fair)
= odds (fair + poor) / odds ( excellent + good + average)
= odds (average + fair + poor) / odds ( excellent + good)
= odds (good + average + fair + poor) / odds ( excellent)
=e^(-1.46)=0.23 （更高维修状况等级为reference，在更低维修状况的odds为0.23）
stata实际进行的分析操作 不够直观
此处stata已经帮我们转化成 比较差的维修状况=0，为reference
=e^(1.46)=4.29（更低维修状况等级为reference，在更高维修状况的odds为4.29）
进口车和国产车相比，在“更低的车辆维修状况等级”的odds 是在“更高维修状况等级“的0.23倍
进口车和国产车相比，在“更高的车辆维修状况等级”的odds 是在“更低维修状况等级“的4.29倍
ologit rep77 foreign, or
一定要记住stata给出的OR是以较低的为reference
ologit rep77 foreign length mpg, or
在控制了汽车的长度、里程之后，进口车有着更高车辆维修状况等级的odds是国产车的18.12倍（95%CI：3.85，85.32）
在控制了汽车的产地、里程之后，车辆每增加1 inch，有更高车辆维修状况的odds增加8.64倍（95%CI：3.90，13.58）
mpg？如何解释
窗口菜单
statistics → ordinal outcomes → ordered logistic → regression
第二十五讲 无序多分类Logistic回归 logistic regression(multinomial)
sysuse fullauto, clear
ologit rep77 foreign, or
进口车（foreign=1）有着更高车辆维修状况等级 的odds是国产车（foreign=0）的4.29倍（95%CI：1.51，12.13）
安装gologit2命令
gologit2命令
- 满足proportional odds假定
gologit2 y x1 x2 x3…, pl or(pl 又叫parallel line 假定)
这个command 和 ologit command 给出的结果相同
- 不满足proportional odds 假定
gologit2 y x1 x2 x3…, npl or
- 检验是否满足proportional odds假定
likelihood-ratio test: lrtest(检验p值是否小于0.05)
gologit2 rep77 foreign, pl or
原理：5个level分成4个二分类变量
category>poor vs cat≤poor
cat>fair vs cat ≤fair
cat>average vs cat ≤average
cat>good vs cat ≤good
有序多分类 当proportional odds假定成立时
进口车（foreign=1）和国产车（foreign=0）比：
odds (excellent + good + average + fair) / odds (poor)
= odds ( excellent + good + average) / odds (fair + poor)
= odds ( excellent + good) / odds (average + fair + poor)
= odds ( excellent) / odds (good + average + fair + poor) =4.29
gologit2 rep77 foreign, npl or
3.94*10^7这个estimate是很不准，上节课 有一个level-excellent的车辆数量是0的，所以estimate很大，standard error也很大，但并不影响interpret results
同样是 category>poor vs cat≤poor
cat>fair vs cat ≤fair
cat>average vs cat ≤average
cat>good vs cat ≤good
有序多分类变量 当proportional odds假定不成立时
进口车（foreign=1）和国产车（foreign=0）比：
odds (excellent + good + average + fair) / odds (poor)=0.93
odds ( excellent + good + average) / odds (fair + poor)=3.45
odds ( excellent + good) / odds (average + fair + poor)=3.28
odds ( excellent) / odds (good + average + fair + poor) =3.94*10^7
检查proportional odds假定是否成立
H0：non-proportional odds 模型可以更好解释结局变量各个等级之间关系
gologit2 rep77 foreign, pl or
est store A (将上面的model存到A中
gologit2 rep77 foreign, npl or
est store B
lrtest A B
p>0.05，拒绝H0，non-proportional odds并没有更好解释结局变量各个等级之间关系，因此可以说满足proportional odds假定
以上为有序多分类变量不满足proportional odds假定的情况，接下来为无序多分类logistic回归
无序多分类logistic回归
把结局变量的某个分类作为reference，然后比较结局变量其他分类相对于reference的相对风险（relative risk）
j指的是结局分类其他变量
有序和无序多分类比较
有序多分类logistic回归：ORj =Pr (cat > j) / Pr (cat ≤ j)
ologit y x1 x2 x3…, or
无序多分类logistic回归：RRj =Pr (cat = j) / Pr (reference cat)
mlogit y x1 x2 x3…, rrr baseoutcome ( j ) ( baseoutcome 不是必须的，若不指定，stata会自动给一个reference level进行比较
mlogit rep77 foreign, rrr baseoutcome (1) poor=1即baseoutcome
无序多分类
进口车（foreign=1）和国产车（foreign=0）比：
Risk (fair) / Risk (poor) =0.20
Risk (average) / Risk (poor) =0.70
Risk (good) / Risk (poor) =1.08
Risk (excellent) / Risk (poor) =1.32*10ˆ7
窗口菜单
statistics > categorical outcomes > multinomial logistic regression
注意gologit2非stata自带command 因此没有窗口菜单操作
26-28讲：数据清理（包括数据导入
29-30讲：高级绘图
想学的没有学到？
第二十六讲 数据整理1 intro to data management 1
26-28讲，我们将使用Do-File（第十讲 第十一讲）
Stata的须知
Stata区分大小写！Stata is case-sensitive（无论是变量的名字 还是command
常见符号
= 是 赋值（e.g. gen x=1
== 是 恒等（e.g. gen x=1 if y==2
| 是 或者（e.g. gen x=1 if y==2 | y==3
& 是 且（e.g. gen x=1 if y==2 & z==3
数据整理常见步骤（stata command）
查看工作路径（pwd）；改变工作路径（cd）
导入excel文件（import excel），CSV文件（import delimited），dta文件（use）
添加变量或变量数值标签（label）
生成新变量（gen）或统计量变量（egen）
讲观测值按照变量数值大小排序（sort；gsort）
改变变量前后顺序（order）
将数据集进行长宽转换（reshape）
合并数据集（merge）
删除（drop）或保留（keep）观测值或变量
导出excel，CSV文件（export）或dta文件（save）
复习stata中添加注释的方法：stata中绿色的文字都是不运行的
三种添加注释的方法，见do file 26
do clear all 注意 // 与all之间需要有空格，否则会报错
Stata IC最多只有2,048个变量; SE:32,767; MP: 120,000
尽管是MP，但一开始打开stata时，默认5000给变量，因此：如果你的变量数>5000，需要设置。例如: set maxvar 8000
正式数据分析开始之前，我们可以打开一个log文件，这样可以将屏幕上的东西都保留在log文件里，便于与co-author其他作者进行分享
若没有指定，自动保存在工作路径下
use "child.dta", clear 该command没有指定去哪儿找child.dta 因此会去pwd工作路径下寻找
use child.dta, clear 双引号可以不加，但如果文件名中有空格，就必须加“”双引号，因此建议都加双引号“”
import delimited "child.csv", clear 报错了！！！ 以前用的方法 在Index of /RePEc/bocode/i (bc.edu) 下载import command 也没能解决 因此有了下面的链接的方法
stata17运行遇到“Java installation not found” - 知乎 (zhihu.com)
安装好stata17后遇到Java installation not found问题解决方案，亲测可行
1、下载java17以上的版本：https://www.oracle.com/java/technologies/javase/jdk17-archive-downloads.html
（注意：下载压缩包版本，windows64位的不建议下载jbk.17.0.6版本，因为很慢）
2、将压缩包解压后将文件放到stata17安装的路径下
3、在stata中输入命令（路径地址根据自己java保存的地址更改）
set java_home "C:\Users\ramfa\Desktop\jdk-17_windows-x64_bin\jdk-17.0.2\"
4、指令运行结束后，关闭stata，再打开就好了
方法来源于简书作者陆陆柒柒https://www.jianshu.com/p/95fe650fb
上述文件已下载在：D:\stata17\外部命令\jdk-17.0.6_windows-x64_bin
import delimited "child.csv", clear
此处stata很聪明的把变量名导入了，没有把它当作第一行的observation，可以看到变量名都是小写的
在导入CSV文件时，Stata自己判断是否把CSV文件的第一行作为变量名导入（但是有时会出错 因此更保险的方法如下：指定stata不把csv的第一行作为变量名导入 import delimited "child.csv", varnames(nonames) clear
我们也可以让Stata必须把CSV文件的第一行作为变量名导入 import delimited "child.csv", varnames(1) clear
导入delimited file时，stata会自动把变量名都变成小写；我们可以让Stata导入的时候把变量名全变成大写 import delimited "child.csv", varnames(1) case(upper) clear
其他内容见Lecture 26.do 文件
以上代码很多 记不住怎么办？可以窗口菜单操作
窗口菜单 运行的命令，如果想要保存至log中，下次直接用：import delimited "D:\stata17\外部命令\数据导入-stata教程-医咖会\child.csv", delimiter(comma) varnames(1) case(upper) clear
log.file中一行太长不好看，换行的话，可以用 ///，记得///和command之间需要空格
在导入Excel文件时，Stata默认第一行不是变量名…余下见lecture26.do文件
import excel "child.xlsx", clear
import excel "child.xlsx", firstrow clear
也可以使用窗口菜单导入excel文件
其他形式数据怎么导入？窗口菜单（本节课学了dta csv excel）
还有其他不支持的格式 可以用 stat/transfer 或者 R语言进行数据转换，一般转化成csv文件是最稳妥的，因为不管什么统计学软件都可以读取csv
第二十七讲 数据整理2 intro to data management 2
依旧使用lecture26
import excel "child.xlsx", firstrow case(preserve) clear
以CHILD_SEX变量为例，进行数据的清洗和整理
codebook CHILD_SEX 字符型变量 因为“1” “2” “M”有双引号
option：replace 在destring时 replace原有变量或generate新变量
给变量和变量的数值加标签
注意：都有空格
生成新变量 变成二分类变量
* 特别注意：在Stata中，缺失被定为为无限大
观测值排序
sort 从小到大 gsort + 从小到大 gsort -从大到小
第二十八讲 数据整理3 intro to data management 3
依旧用lecture26.do
打开browser 可以直接点击 也可以browse
变量排序
合并数据集
child.dta 每个孩子一条观测值，也就是每个孩子只有一行，现在我们想合并anthro.dta，为一个孩子的各种数据的dta，每个孩子有好几条观测值，好几行数据，分别是不同的时间点进行的测量，是一个长的long dataset
打开一个数据集，打开的就叫master dataset
master dataset（child） using dataset（anthro）
one-to-one 一行 一行
many-to-one 多行 一行
one-to-many 一行 多行（本次就是第三种情况）
many-to-many 多行（未知） 多行（未知）（自己不知道多少行，让stata自己决定，耗时长，并且不佳，知道的话尽量还是对应用以上三种之一）
3580未merge，并且来自master，说明anthro是child的子集，而child是初始baseline，因此有孩子失访了
自动生成了一个新变量_merge，master only即只在master中，不在using dataset中
只想保留match成功的（既有maternal information又有infant anthropology measurement的孩子）：keep if _merge == 3 drop _merge 完成数据集合并
reshape 有时我们有一个长数据集，想变成宽数据集，反之亦然
help reshape 理解：i是id，j是时间点，stub为测量内容（如身高），长数据集：时间点是分开的，宽数据集：时间点跟在stub后面
child数据集中 j对应visit（各个时间点）
reshape long stub, i ( i ) j ( j ) （reshape (to) long变成长数据集）
reshape wide CHILD_ADJ_AGE WEIGHT-HEAD_CIRC, i(CHILD_PIDX) j(VISIT)
-意为从x1变量到xn变量
long → wide 变量名增加了，visit被drop掉了
8096为2024的四倍
因为：每个孩子的visit不一定都有四行，long→ wide→long，还原的时候stata自己补充了
egen（egenerate）生成统计量变量
by CHILD_PIDX: egen max_visit = max(VISIT)
egen height_mean = rowmean(HEIGHT6 HEIGHT12 HEIGHT24 HEIGHT36) 四个时间点的身高平均值
egen height_miss = rowmiss(HEIGHT6 HEIGHT12 HEIGHT24 HEIGHT36) 看这四个变量中，每个孩子到底有几个缺失值
有没有确实少一点的？tab height_miss
egen height_nonmiss = rownonmiss(HEIGHT6 HEIGHT12 HEIGHT24 HEIGHT36) 看一下没缺失的有几个 height_miss 和 height_nonmiss 互补
help egen 自己探索！很强大的功能 多看看下面的examples
export excel "Sample", replace //有问题！
没有保存变量名字（第一行）！
export excel "Sample.xlsx", firstrow(variables) replace //“变量名”作为第一行导出
excel打开csv文件，乱码，但其实数据是完好的，因此最好把变量和变量标签都保存成数值格式，界面也尽量使用英文的
窗口菜单操作
可以勾选是否保留第一行名字等
第二十九讲 高级绘图1 advanced graphing 1
第8行 保留变量keep 保留观测值keep if ；保留了所有matched观测值child
9… 去掉 变量_merge 多余的行，因为这一行的作用就是筛选下_merge==3的child
13…
14…
15…不知道蓝/红色是谁，待会教加标签
22…
换不起行？因为///前没有空格，stata把///当成command了
√
缺点：都排在一列上了，想让他们更分散均匀点
option：jitter ( ) 抖点-让点散开
twoway (scatter WEIGHT MOM_AGE if CHILD_SEX==1, jitter(2)) ///
(scatter WEIGHT MOM_AGE if CHILD_SEX==2, jitter(2))
缺点：单个点太大了，颜色可以接受
twoway (scatter WEIGHT MOM_AGE if CHILD_SEX==1, jitter(2) msize(tiny)) ///
(scatter WEIGHT MOM_AGE if CHILD_SEX==2, jitter(2) msize(tiny))
twoway (scatter WEIGHT MOM_AGE if CHILD_SEX==1, jitter(2) msize(tiny)) ///
(scatter WEIGHT MOM_AGE if CHILD_SEX==2, jitter(2) msize(tiny)) ///
(lfit WEIGHT MOM_AGE if CHILD_SEX==1) ///
(lfit WEIGHT MOM_AGE if CHILD_SEX==2)
缺点：颜色有黄有绿，更不知道谁是谁了
缺点：置信区间有点细，想加粗 clwidth(thick)
观察到x轴在15之后，45之前出现点，想改一改x轴
xlabel(15(5)45)
y轴想延长到25 ylabel(5(5)25)
xtick(15(1)45) ytick(5(1)25) 间隔1 画tick，就是坐标轴的小尺度
做好一个图后，更改细节就很快很方便了
第三十讲 高级绘图2 advanced graphing 2
学习一个新的command：coefplot
如果之后的文章使用了coefplot，最好引用一下这篇文章
ssc install coefplot 报错了！
cannot write in directory C:\Users\���\ado\plus\c
解决办法：https://bbs.pinggu.org/thread-10685955-1-1.html
每次记得运行一下这个do.file 就可以使用ssc install 命令了
第9行
图上给出了这几个的β estimate 以及95% CI
但在实际中，我们很少关注截距是怎样的，而且截距的范围很大（-5000，20000）。导致mpg trunk等置信区间被压缩成一个点
第13行 coefplot, drop(_cons)
第15行 coefplot, keep(trunk turn)
第32行 coefplot D F, drop(_cons) xline(0) 缺点：legend不清楚D、F代表什么 给它们加label
灰线和黑线均匀分布在纵坐标对应刻度的两边，两条线的距离可以更改
1万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
