Part1

第一页致敬素未谋面的无冕之王朋友，原因在这！“关于“真实世界数据公众号MSTATA软件是否涉嫌版权侵权的探讨.

小编初心是身边很多医生朋友都苦于好的rcs图和讨论，于是好心分享给大家免费用于学术的rcssci包，这个包是经过自己多篇文章和心血总结开发，小编发布rcssci的时候我就特地选择了商用授权的license(Artistic2.0)，以避免后期可能发生纠纷。但是居然被人潜在商用，且未经任何本人许可，违背该包商用原则，及全部抄袭公众号原创思维未加任何引用说明。

请问你在全网看到别的rcs系列原创有实战统计的思维内容吗？别说R包调用问题，人家rms ggplot2开放的MIT自由协议允许你商用的你随便改，rcssci这个是不允许商用需要授权，商用你懂版权吗？可以说有些人毫无版权意识，视协议于无物，加上后台公众号主留言其涉及多处完完全全抄袭更加实锤其行为。别人辛苦开放给学术，不是让你拿去未经许可商用的！

留言中是否有其相关人士转移话题洗白，我建议其大可不必，群众眼睛是雪亮的，反正至今仍然无人联系我授权问题，如果涉及商用小编会继续跟进相应事宜走诉讼程序让对方知道什么是侵权。当然，你免费开放小编十分欢迎，毕竟造福广大医生群体，但是依然需要原作者授权引用出处来源！

鉴于此恶劣事情，小编十分寒心，决定以后都不会再公开发布rcssci包及完成的其他便捷包核心代码，以免此类事件再发生。

Part2

下面开始今天正文。

本文内容很多很硬，需要阅读前面rcs系列文章、人年发病率文章作为基础铺垫。这是一篇迟到8年的总结，正好近期写论文遇到几篇队列研究的顶刊都发了纵坐标为incidence rate (per 100 patient-years at risk)的限制性立方样条RCS，痛苦的小编不得不打起精神给自己方法学版图填上8年来的空白。

一、HR、RR、incidence rate的y轴取舍

心血管或肿瘤等领域，研究设计常见的队列cohort，随访时间(time)通常较长，且每一位患者time可能不同。研究者不仅仅关心结局发生与否，更关心结局什么时间发生，即干预方法能多大程度上延长患者生存时间。因此，这一类结局指标，又叫做时间事件（time to event）变量，event 即结局发生与否，time即结局发生的时间(阳性)或随访全周期时间(阴性)。在做回归分析时候，根据time是否精准或者能否获取，可以采用poisson或quasipoisson进行RR值、cox分析进行HR值等不同效应量估计。随着方法学的进步，RCS孕育而生。尤其是剂量效应非线性关系，从简单线性的x-y到广义线性的g(x)-y到孟德尔snp-y，都有RCS的影子。本文局限在time相关队列的RCS，考虑x与y相对率指标的效应量推荐等级是HR = RR > OR，绝对率指标的效应量incidence rate (per-year) > incidence rate。

1.1 HR

范文1(PMID: 36094529)，采用的是cox回归RCS，y轴为标准的HR值。原文中在每天步行2000步的时候，全因死亡风险对数值(lnHR)降低8%，随着步数的增加，风险持续降低，直到每天10000步。此外，心血管疾病和癌症发病率也有类似的关联。https://mp.weixin.qq.com/s/DSQ-SbKy0MlfGPb1ujBNlw 文章解释有误，需要注意。本范文RCS图也有误，y轴坐标应该标为lnHR。

1.2RR

时间关联的队列剂量效应分析，rcs采用RR指标。rcs描述：To examine the dose-response relationships between beverage intake and the outcomes, we used restricted cubic spline regression with three knots. Tests for non-linearity were based on the likelihood ratio test comparing two models: one with only the linear term and the other with the linear and the cubic spline terms. 最终rcs的y轴是RR而不是HR，可能是因为ph假设不通过，cox分析采用的是时间依赖性cox(Time varying Cox proportional hazards models, conditioned on age and follow-up cycle, were applied to estimate hazard ratios and 95% confidence intervals for the associations of each beverage intake with all cause mortality, )，作者退而求其次就做了RR的RCS？

1.3人年发病率

小编团队近期参考的两篇顶刊范文，均采用RCS y轴为incidence rate( per100 patient-year)。那么，此时你是否有困惑(至少小编起初是有的)，为什么有time了不用HR做y轴呢? 因为发病率和RR回答的角度不同！RR回答相对风险倍数，incidence rate回答绝对发病率。

RR说明的是暴露（药物）与结局（死亡）的关联强度及方向。关联强度是指试验组研究对象的死亡率是对照组研究对象的“多少倍”，RR值越偏离1，说明“倍数”越大；而方向，是指正关联（RR>1）还是负关联（RR<1）。比如RR=0.5，代表药物和心脏病患者死亡之间存在负向的关联，A药和B药均可以将死亡率降低到原来的0.5倍（即下降50%）。对于A药和B药，两个试验的RR都是0.5。我们可以说“A药可以使室性心律失常患者的死亡率下降50%；B药也可以使急性心肌梗死患者的死亡率下降50%”，相对危险度降低（Relative Risk Reduction, RRR），RRR=1-RR。但RR和RRR不能说明的是，相比于对照组，药物/干预方式降低试验组研究对象死亡率的绝对值是多少，即可以让多少比例的患者免于死亡。

二、人年发病率poisson

2.1 poisson原理

Poisson回归一般用于单位时间、单位面积、单位空间内某事件发生数的影响因素的探讨，事件的发生服从Poisson 分布。如果某一事件的的观测频数服从一个均值为λ的Poisson分布，则Poisson回归模型如下：：λ=e^Xβ，ln?λ=β_0+β_1 x_1+β_2 x_2+?β_m x_m β0=ref表示各个“自变量取值为0”时观测频数的自然对数值，注意这里的取值为0的水平指的是ref水平，spss为最高level。βi则表示其他自变量取值不变时，自变量xi每改变一个单位，所引起的观测频数λ的自然对数值的改变量。当不同观察单位的发生事件观测基数不一样时，也就是临床长期随访最常见的情况，需要将发生数除以相应单位的总观测数，用于抵消观察单位数的不同对结果造成的的影响。

2.2 二手归纳数据的poisson

Poisson计算很多教程全是二手数据即归纳好的变量组合和人年，小编也被坑了很多年，不懂为什么教科书这么多年永远是一样的案例?其实，原始数据的poisson才是更靠谱的更接地气的，不然折腾半天搞个汇总有多大应用价值？对于二手总结带人年数据，SAS实现且解释好，推荐。SPSS有点蠢(傻变量最高level为ref参考值，反向解释)且在不饱和模型中基线率可能不正确。以下案例来自胡良平老师《SAS统计分析教程》。

2.2.1 二手数据SPSS对数线性模型

2.2.2 二手数据SAS genmod

2.2.3 二手数据R

二手数据R需要glm即可完成，绝大时候仅仅用来计算统计推断的RR95%CI。

2.3 原始数据的poisson

2.3.1 原始数据计算RR或人年发病率原理

这里要明确一个理念，实际人年incidence rate在上一章节阐述过，而做RCS是建立quasipoisson回归生成个体预测人年发病率 (pre.率)及95%CI再拟合校正多因素的RCS曲线。

2.3.2 SAS genmod

见rcssci0.4.0推文中介绍，采用“队列研究RR值.sav”数据进行SAS9.4与R4.2.的结果一一比较。RR标准计算是广义线性模型的log函数link，有公众号提供了phreg命令但是结果错误(神坑！)。Genmod的SAS 结果稳定可靠，可见案例文件hp的RR=1.77(1.01-3.12)。而若是采用logistic分析，OR=2.44(1.05,5.70)，置信区间更宽，点估计值更大了。R结果采用quasipoisson结果与SAS一致，而采用poisson会结果偏倚。笔者检索后发现，拟合模型的残差比自由度大很多，这个时候要使用quasipossion，即quasipossion更稳定， RR=1.77(1.04-3.15)与SAS结果几乎一致。

2.3.3 R实现个体人年发病率RCS

glm和rms::Glm函数都可以实现人年发病率预测模型预测概率，后者可以无缝连接rcs。rms包极其强大，内置了8种类可以制作RCS且支持Predict函数的家系模型，有ols、lrm、cph、psm、bj、Glm、Gls、Rq，琳琅满目。但是细节也很多，稍微不慎就容易踩坑，比如上一期介绍的svycoxph。

前面案例列举了基于R的二手数据计算RR及可预测分层人年发病率，但是使用场景很罕见，且基于二手汇总的poisson预测效果不佳。所以，实际情况我们都是基于个体的数据预测建模，更加实用预测效果也更好。尤其是x变量是连续型例如age，相对于年龄分类型agec变量，age可以很好的制作平滑的人年发病率RCS！注意避坑坑多的可怕如履薄冰，具体细节不再公开提供解决方案！嗯，下面就是人年发病率rcs闪亮登场的时刻！

小结

“

队列的时间time记录细节不同方法学不同，请参阅思维导图。 

人年发病率RCS关注绝对率侧重描述，RR/HR RCS关注相对率(风险倍数)侧重效应推断。

”