复杂性状的细胞类型表达特异性整合（CELLECT）是一个计算工具包，用于识别复杂性状潜在的可能病因细胞类型。CELLECT利用现有的遗传优先级模型，在识别可能的病因细胞类型时，整合单细胞转录组和人类遗传数据。
工作原理
CELLECT使用已建立的遗传优先级模型，如S-LDSC和MAGMA协变量分析，量化了常见的多基因GWAS信号（遗传力）和基因的细胞类型表达特异性（ES）之间的关联。CELLECT的输出是给定人类复杂疾病或特征的优先病因细胞类型列表。
CELLECT将GWAS数据和细胞类型表达特异性估计作为输入。为了计算ES的稳健估计，我们开发了一种称为CELLEX（CELL型表达特异性）的计算方法。CELLEX建立在不同ES指标提供互补的细胞类型表达特异性图谱的观察基础上。我们的方法结合了“群体智慧”方法，通过整合多个ES指标来获得改进的稳健性和更具表达性的ES指标，该指标可以捕捉表达特异性的多个方面。
安装（centos7 env）
Step 1: Install git lfs
sudo yum install epel-release)
sudo yum install git
curl -s https://packagecloud.io/install/repositories/github/git-lfs/script.rpm.sh | sudo bash
sudo yum install git-lfs
git lfs install
Step 2: Clone CELLECT repository
git clone --recurse-submodules https://github.com/perslab/CELLECT.git
Step 3: Install Snakemake via conda
conda install -c bioconda -c conda-forge snakemake">=5.27.4"
conda install -c conda-forge mamba
mamba create -c conda-forge -c bioconda -n snakemake snakemake
配置文件
CELLECT 是利用snakemake进行一键化运行的程序，config.yml是CELLECT提供给我们的一个snakemake配置文件，该文件中规定了我们需要的输入文件，输出路径，特殊关注的细胞类型等信息，改文件中有详细的注释，可以自行阅读。
输入文件1：SPECIFICITY_INPUT:
该文件是基因的细胞类型表达特异性文件，CELLECT给出了两个实例文件，example/mousebrain-test.csv；example/tabula_muris-test.csv，针对我们自己的scRNA数据，可以使用CELLex工具得到该文件。具体代码如下：
import numpy as np
import pandas as pd
import cellex
data = pd.read_csv("./data.csv", index_col=0)
metadata = pd.read_csv("./metadata.csv", index_col=0)
eso = cellex.ESObject(data=data, annotation=metadata, verbose=True)
eso.compute(verbose=True)
eso.results["esmu"].to_csv("mydataset.esmu.csv.gz")
其中我们得到的"mydataset.esmu.csv.gz" 文件，即分析需要的输入文件。注意，如果我们的单细胞数据基因使用的是symbol_name，需要替换为ensem_gene_ID。同时需要注意的是，此文件中celltype名，中间不能有空格，需要用下划线进行连接。
输入文件2：GWAS_SUMSTATS:
该文件是多基因GWAS信号，需要我们根据自己的分析提供，如果没有相关分析，可以再GWAS数据库（https://gwas.mrcieu.ac.uk/）中下载疾病相关的vcf文件。
在拿到最终的输入文件之前，CELLECT为我们提供了一个整理脚本，实例文件下载方法：
wget https://portals.broadinstitute.org/collaboration/giant/images/c/c8/Meta-analysis_Locke_et_al%2BUKBiobank_2018_UPDATED.txt.gz -P example/
#--- 格式
CHR POS SNP Tested_Allele Other_Allele Freq_Tested_Allele_in_HRS BETA SE P
N
7 92383888 rs10 A C 0.06431 0.0013 0.0042 7.5e-01 598895
12 126890980 rs1000000 A G 0.2219 0.0001 0.0021 9.6e-01 689928
4 21618674 rs10000010 T C 0.5086 -0.0001 0.0016 9.4e-01 785319
#---
# 可以查看该文件的格式将我们自己下载的GWAS数据进行整理，
python ldsc/mtag_munge.py \
--sumstats example/Meta-analysis_Locke_et_al+UKBiobank_2018_UPDATED.txt.gz \
--a1 Tested_Allele \
--a2 Other_Allele \
--keep-pval \
--p PVAL \
--merge-alleles data/ldsc/w_hm3.snplist \
--out example/BMI_Yengo2018
# ---------------- or ----------------------
wget https://www.dropbox.com/s/ho58e9jmytmpaf8/GWAS_EA_excl23andMe.txt -P example/
# ---
MarkerName CHR POS A1 A2 EAF Beta SE Pval
rs13090388 3 49391082 C T 0.6905 -0.02852 0.00184 4.29e-54
rs7630869 3 49522543 C T 0.6922 -0.02848 0.00184 4.61e-54
# ---
python ldsc/mtag_munge.py \
--sumstats example/GWAS_EA_excl23andMe.txt \
--merge-alleles data/ldsc/w_hm3.snplist \
--n-value 766345 \
--keep-pval \
--p PVAL \
--out example/EA3_Lee2018
如上输出的BMI_Yengo2018开头的文件即为我们需要的另外一个输入文件。此处注意，vcf文件有些提供的是LP（-log10 pvalue） 值，需要转化到Pvalue值；第二，如果我们在文件中没有找到需要的统计量，可以进行推导计算，具体计算方法如下：
# 1. se 计算
## 1.1 有 beta、p ，计算se；
se=sqrt(((beta)^2)/qchisq(p,1,lower.tail=F))
## 1.2 有 or、p ，计算se；
se=abs(log(or)/qnorm(p/2))
# 2. zscore计算
## 2.1 有 beta、p ，计算zscore；
c=qnorm(1-p/2)
### 情况1：beta > 0，zscore = c
### 情况2：beta < 0，zscore = -c
## 2.2 有 or、p ，计算zscore；
c=qnorm(1-p/2)
### 情况1：or > 1，zscore = c
### 情况2：or < 1，zscore = -c
## 2.3 有 beta、se ，计算zscore；
zscore=beta/se
## 2.4 有 or、se ，计算zscore；
zscore=log(or)/se
程序运行
文件准备好，填写路径到config.yml文件中，就可以进行程序运行了，两种模式的运行命令如下：如果以上提到的坑大家没有掉下去，此处一般也不会出现问题。
# CELLECT-LDSC:
snakemake --use-conda -j -s cellect-ldsc.snakefile --configfile config.yml
# CELLECT-MAGMA:
snakemake --use-conda -j -s cellect-magma.snakefile --configfile config.yml
最后我们在指定的路径下成功获取到了我们的分析结果：results/prioritization.csv，需要绘图的也可以根据我们的结果进行绘图了。
gwas specificity_id annotation beta beta_se pvalue
bbj3 EAT_set1 T_cell 2.4251659065799996e-10 1.6731792749e-10 0.07360796137730001
bbj3 EAT_set1 Endothelial 1.65815028102e-10 1.16935618762e-10 0.0780949820402
bbj3 EAT_set1 ADIP 1.9396405056900004e-10 2.0653247161299997e-10 0.17382801157699998
