前言

资深杰青韩敬东团队在前段时间发表了Inference of differentiation time for single cell transcriptomes using cell population reference data，这是一篇关于利用plsc（偏最小二乘主成分回归）来推断细胞发育时间的文章

偏最小二乘回归

基于R语言的偏最小二乘回归可以参考文献《偏最小二乘建模在 R 软件中的实现及实证分析》
那么我们知道最小二乘法往往针对于一个因变量对多个自变量的回归；二偏最小二乘通常针对的是多个因变量对多个自变量的回归

偏最小二乘回归

每一个因变量对应一个线性模型

那么偏最小二乘回归的步骤为：

1.将因变量矩阵（Y）和自变量矩阵（X）进行数据中心化并进行PCA分解，对于因变量矩阵（Y）分解的主成分为：u1，u2，...，up；而对于自变量矩阵（X）分解的主成分为v1，v2，...，vq

分解如下图所示：

PLS的核心思想：

2.利用典型相关分析的思想，寻找u与v中相关性最大且贡献也最大的一对PC向量（u1和v1），用u1和v1分别做因变量矩阵（Y）和自变量矩阵（X）的线性回归，若回归方程达到响应的精度，则停止算法；若没有达到精度，则利用余下的信息在进行一轮分解（这里的余下信息理解为 Y - w * u1 和 X - c * v1 所得的矩阵再次进行分解，再次寻找相关性最大的一对PC向量计算回归方程，最后将每一次所求的回归方程进行加和），最终还原出 y1...yi 与 x1...xj 的回归关系

pls与单细胞轨迹推断

作者的测序方式是对不同时期不同种类（每种种类的细胞有生物学重复）的单细胞做了"bulk-seq"，并以基因作为自变量的指标

cpdata
而因变量为：
cptime
事实上因变量只有一列数据，表示不同细胞类型的发育（取样）时间
脚本下载自：https://www.picb.ac.cn/hanlab/iCpSc.html

library(pls)

cpdata = read.table("CellPopulation.DEGs.log10_overlap.txt",sep='\t',header=T,row.names=1)
cptimes = as.numeric(read.table("CellPopulation.DEGs.log10_overlap.txt",sep='\t',header=T,row.names=1)[1,])

t_cpdata = as.matrix(t(cpdata))
t_cpdata_frame = data.frame(cptimes, t_cpdata=I(t_cpdata))
# pls建模
mod = plsr(cptimes~t_cpdata,data=t_cpdata_frame,ncomp=6, scale = TRUE, validation = "none")

对于X矩阵与cptimes（Y 矩阵）来说， 1 comps 为X矩阵经过PCA分解后的PC向量中（v1，...，vq）与 Y矩阵PCA分解后的PC向量中（u1，...up）相关性最强且贡献最大的那一对PC向量（假设u1和v1为相关性最强且贡献最大的那一对PC向量，典型相关分析的思想），其余的 comps 是利用余下的信息PCA分解后重复上述分析所得

余下信息解释：

因此plsc分解后可视化出来的其实是带有一定的基因表达信息：

由于每个时间点细胞中基因表达不同，因此可以很好的分开，从而显示这种细胞间发育时间的关系（PLSC1，PLSC2即为上述的 1 comps，2 comps ），作者只利用了1 comps，2 comps来做回归，来还原y1...yi 与 x1...xj的回归

并根据基因表达特征预测发育时间：

predT=predict(mod,ncomp=2,type='response',newdata=t_cpdata)

总结

为了避免共线性，作者利用PLS建立了基因每个样品的基因表达特征与每个发育时间的线性关系，以此来根据基因的表达特征来预测所处于的发育时间