童筱君，李辰璐，王國芬，朱小春，王曉冰

1.臺州市中心醫(yī)院（臺州學院附屬醫(yī)院） 風濕免疫科，浙江 臺州 318000；2.溫州醫(yī)科大學附屬第一醫(yī)院 風濕免疫科，浙江 溫州 325015

強直性脊柱炎（ankylosing spondylitis，AS）是一種慢性炎癥性自身免疫性疾病，常導致脊柱疼痛、僵硬及終末期關節(jié)融合。AS的發(fā)病率依種族而異，在中國約為0.3%[1]。AS病因未明，被認為與遺傳和環(huán)境相關[2-3]，尤其是遺傳因素在AS發(fā)病中起著重要作用，其遺傳度超過90%[4-5]。

人類白細胞抗原（human leukocyte antigen B27，HLA-B27）是迄今為止所知的復雜疾病中相關性最強的一個等位基因[6]。國際AS遺傳學聯(lián)盟發(fā)現(xiàn)除了HLA-B27外，還有31個不同的基因位點被確定與AS相關，大部分位于主要組織相容性復合體（major histocompatibility complex，MHC）區(qū)域。四種編碼位點的氨基肽酶在MHC I-肽復合物抗原遞呈通路中發(fā)揮作用[7]，可能在AS的發(fā)病機制中起重要作用。國內外已有大量關于ERAP1 和ERAP2 兩個氨基肽酶編碼位點的研究，以及將相關研究結果轉化為臨床治療的應用[7-8]。而另一個與疾病相關的氨基肽酶編碼位點NPEPPS 在AS中的后續(xù)研究尚不足。本項目利用免疫芯片的數據和二代RNA-Seq技術對疾病相關NPEPPS基因的遺傳變異進行eQTL分析，探討NPEPPS在AS發(fā)病和進展中的作用。

1 對象和方法

1.1 對象 選取溫州醫(yī)科大學附屬第一醫(yī)院風濕免疫科收治的54例AS患者作為AS組，78例健康對照者為對照組。AS組所有患者均符合1984年紐約分類標準[9]。本研究經溫州醫(yī)科大學附屬第一醫(yī)院倫理委員會批準，所有患者簽署了知情同意書。

1.2 免疫芯片基因定型 從入選者5 mL全血中分離出外周血單個核細胞（peripheral blood mononuclear cells，PBMCs），提取DNA，用Illumina平臺進行單核苷酸多態(tài)性（single nucleotide polymorphism，SNP）微陣列的基因定型檢測，基因型被定位到人類基因組NCBI Build 37（hg 19）。用已知的基因型預測未知的基因型并對缺失的數據進行補缺（imputation）和SNP的選擇。樣本的基因型用PLINK格式（ped和map），提供患者ID和基因組每個SNP的堿基。用PLINK提取包含NPEPPS的染色體區(qū)域的SNPs（Chr17，hg 19）。產生每個SNP的次要等位基因頻率（minor allele frequencies，MAF），并用MAF＞0.3來進行SNP的剔除。Map文件上的SNP位點被轉換到hg 19，去除標志SNPs。利用參考基因組的SNP數據，根據鄰近的變異位點間的連鎖不平衡（linkage disequilibrium，LD）來推斷未定型的SNPs。在補缺之前用SHAPEIT和IMPUTE2進行缺失基因型的補缺，用info＜0.5 來去除補缺較差的SNPs。用PLINK計算免疫芯片的基因型，并設定窗包括所有與疾病相關SNPs（rs9901869）相關性R2＞0.1，以便進行后續(xù)的分析。用PLINK抽取這些SNPs，將基因型從等位基因對轉換為最小等位基因計數（0、1或2，分別代表主要基因型、雜合基因型和次要基因型）。剩余的MAF＞0.01的SNPs被分成不同的組，每組都與補缺的SNPs密切相關（即組內R2=1），這若干組不同的SNPs將進入后面的分析。

1.3 RNA測序方法 將RNA從PBMCs樣本中提取出來，反轉錄成cDNA，用Illumina TruSeq標準的全RNA樣本試劑盒預處理，在Illumina平臺進行測定。用TopHat version 2.0.6獲得的短read定位到人類基因組NCBI Build 37（hg 19）。對齊的read應用于HTSeq，產生每個基因的read計數，這可以作為經過read多樣性校正后的基因表達的測定。來源一個基因的不同轉錄產物，差異來源于轉錄起始/終止位點，外顯子、內含子及非翻譯區(qū)的包含與否。用Cufflinks suit來匯編read到融合的個體轉錄文件中，獲得參考轉錄組。用RSEM來將read排列到參考轉錄序列中去，來計算同種型的表達豐度。

1.4 eQTL分析方法 計算正態(tài)化和模型驗證：DESeq2是R語言中一種差異基因表達分析的程序包，用于RNA-Seq計數數據的應用。用內部的函數將計數數據正態(tài)化，并觀察其大小和水平層次，判斷序列的多樣性。轉錄本的計數用DESeq2正態(tài)化，并提取轉錄產物。每個同種型的FASTA序列用Cufflinks的gffread程序獲得。用NCBI網站（Madden，2002）的基礎局部比對工具（BLAST）來跟人類基因組比對，提供外顯子不同或者3’或者5’UTR異常的特征信息。GLME模型用于檢驗基因型對轉錄表達的效應。用Wilcoxon符號秩和檢驗來檢測轉錄表達在AS組和對照組間的差別。基因型和轉錄產物比率用線性模型進行驗證，并對疾病狀態(tài)進行校正。

1.5 統(tǒng)計學處理方法 采用R v3.6.1 for Linux軟件包進行統(tǒng)計學分析。2 組間比較用秩和檢驗，用Wald test計算DESeq2的P值。P＜0.05為差異有統(tǒng)計學意義。

2 結果

2.1 SNP rs9901869 不同基因型對應的NPEPPS中TCONS_00206150轉錄本和NPEPPS的表達比較 對2組進行免疫芯片基因分型和RNA測序的NPEPPS表達水平分析，純合等位基因頻率和雜合等位基因頻率，2組間NPEPPS及其TCONS_00206150轉錄本表達的差異無統(tǒng)計學意義（P＞0.05），見圖1。

圖1 SNP rs9901869不同基因型對應的NPEPPS中TCONS_00206150轉錄本和NPEPPS的表達水平

2.2 2組NPEPPS 不同的轉錄本的表達情況 12種NPEPPS 不同的轉錄本的表達密度見圖2；2組不同的轉錄本的表達水平各不相同，其中2組在TCONS_00206150、TCONS_00206153和TCONS_00206137的表達水平間差異有統(tǒng)計學意義（P＜0.05），見圖3。

2.3 SNP rs9901869不同基因型下NPEPSS基因及其轉錄本表達水平比較 SNP rs9901869 不同基因型所對應的NPEPSS基因及其轉錄本的表達水平各不相同；其中部分轉錄本，如EFCAB13、TCONS_00206150、TCONS_00206157、TCONS_00206153、TCONS_00206137和TCONS_00206130 2組的NPEPSS基因及其轉錄本表達差異有統(tǒng)計學意義（P＜0.05），見圖4。

2.4 SNP rs180677251 基因型的NPEPPS 基因及其轉錄本的表達比較 對主導SNP rs9901869 和附近1 MB范圍內的基因的eQTL分析，發(fā)現(xiàn)SNP rs180677251不同的基因型所對應的NPEPPS基因及其轉錄本TCONS_00206137的表達水平差異有統(tǒng)計學意義（P＜0.05），見圖5。

3 討論

圖2 12種NPEPPS 不同轉錄本的表達密度圖

圖3 2組NPEPPS 不同轉錄本的表達水平比較

到目前為止，尚未發(fā)現(xiàn)與健康對照組相比，AS患者疾病相關的氨基肽酶表達有顯著變化的證據[10-11]。而NPEPPS 是編碼嘌呤霉素敏感性氨基肽酶的基因。NPEPPS 在大多數組織中以胞質氨基肽酶的形式表達，但在大腦中也發(fā)現(xiàn)了一種與膜相關的形式[12-13]，是目前已知的唯一能剪切polyQ序列的胞漿氨基肽酶[14]。鼠模型研究發(fā)現(xiàn)神經退行性疾病中該基因的過表達能使蛋白聚集減少，增強自噬作用[15]。而在AS患者的回腸活檢中發(fā)現(xiàn)自噬的證據，并被證實是IL-23 在腸道表達的驅使者[16]。NPEPPS是否通過自噬作用或者其他機制作用于細胞內肽的處理，仍需要深入探討。至今為止，AS基因表達研究尚未發(fā)現(xiàn)NPEPPS存在基因與基因之間的相互作用[17]。雖然已經發(fā)現(xiàn)一個與AS疾病相關的位于NPEPPS基因的關鍵SNP rs9901869[18]，但是這個SNP對NPEPPS的基因表達是否具有調控作用，以及是否存在其他與AS疾病相關的SNP能影響NPEPPS基因表達，尚不明確。故本項目通過免疫芯片基因分型技術和RNA測序技術，探討AS患者和健康對照者NPEPPS基因上的SNP rs9901869的基因型與NPEPPS基因表達水平之間的關系，發(fā)現(xiàn)AS疾病風險等位基因的基因型與NPEPPS基因及其不同轉錄本的表達水平相關。進一步揭示了NPEPPS的基因位點多態(tài)性影響了NPEPPS的表達水平。鑒于NPEPPS在疾病發(fā)病中的潛在作用，推測NPEPPS的eQTL調控可能與AS的發(fā)病相關。

圖4 不同SNP rs9901869 的基因型下NPEPSS 基因及其轉錄本的表達比較

圖5 SNP rs180677251 基因型的NPEPPS 基因及其轉錄本的表達比較

將基因型信息和基因表達數據進行配對分析，有助于闡明遺傳變異對一個基因轉錄模式的影響，即表達數量性狀基因座。當這些變異也能顯示出疾病相關性時，它們將有助于證明單個SNP在疾病發(fā)病中的作用。如果這些SNPs干擾了那些被細胞器識別并參與mRNA的產生和DNA基序的加工過程，將會影響附近基因的表達（cis-eQTLs），或其他相互關聯(lián)通路上的基因表達（trans-eQTLs）[19-20]。隨著二代測序技術的進展，RNA測序（RNA-Seq）為eQTL分析提供了一個絕好的平臺，使高通量的基因組至轉錄檢測都得以實現(xiàn)。跟傳統(tǒng)的基因芯片不同，RNAseq可以提供更多單基因源轉錄產物的編碼序列，而這個差異可能是由于剪接位點不同、外顯子的包含與否或3’和5’非翻譯區(qū)的交替所引起[21]。這些信息可以檢測到影響mRNA剪接前的那些變異（sQTLs），而這些變異可能造成同一基因在細胞中表達出完全不同功能的異構體。本研究發(fā)現(xiàn)在AS和健康對照者中NPEPSS不同的轉錄本的表達水平各不相同，其中TCONS_00206150、TCONS_00206153和TCONS_00206137的表達水平差異顯著。SNP rs9901869 及其附近的SNP rs180677251的不同基因型所對應的NPEPSS基因及其轉錄本的表達水平各不相同，進一步證實AS中基因型可能影響基因表達的推論。

本研究的缺陷是進行基因研究的樣本數相對較少，雖然這些樣本能夠檢測出大量SNPs，但是仍需要大量樣本的研究證實本研究結果。本研究中發(fā)現(xiàn)AS和健康對照者中AS相關SNP rs9901869與NPEPPS表達的直接關系仍無法確定，仍需要進一步通過細胞和（或）動物模型進一步驗證。

通過對AS和正常人免疫基因的研究，發(fā)現(xiàn)除了HLA-B27以外可能還有其他更復雜的遺傳學背景和環(huán)境因素參與疾病的發(fā)生和進展，比如說與免疫炎癥性疾病相關的SNP有關，以及還有一些與MHC I-肽復合物抗原遞呈通路相關的基因位點變異有關，進一步揭示了AS發(fā)病中參與的重要基因位點和關鍵信號通路，有助于探討AS發(fā)病的潛在機制，為今后的靶向治療提供了依據。