復制性衰老中全基因組染色質開放程度的變化,醫學遺傳學論文摘要：為研究復制性衰老中表觀遺傳學的變化,以年輕代齡(PD26)和年老代齡(PD55)的人胚肺二倍體成纖維細胞(2BS)為材料,進行ATAC-seq測序,并使用Bowtie2、FastQC、MACS、deepTools、phastCons、R語言、HOMER以及基因本體論(GeneOntology,GO)對測序結果進行深切進入分析和挖掘。結果顯示,2BS細胞染色質開放區域的保守性較強,其主要集中于啟動子(promoter)區至轉錄起始位點(transcriptionstartsite,TSS),且年輕細胞與年老細胞的啟動子區至轉錄起始位點的開放程度不同,細胞隨著衰老在啟動子區至轉錄起始位點的開放程度逐步下降。進一步對年輕和年老細胞之間存在顯著差異的靶基因進行GO分析發現,這些差異基因主要與細胞進程、代謝、生物性調節、結合功能、催化功能相關。本次研究發現衰老經過中染色質開放區減少、轉錄調控水平下降,提示復制性衰老與轉錄調控關系密切。本文關鍵詞語：復制性衰老;ATAC-seq測序;轉錄調控;染色質開放程度;Abstract：Tostudytheepigeneticchangesinreplicativesenescence,theATAC-seqwasperformedinyoung(PD26)andsenescent(PD55)2BScells.Aftersequencing,theresultswereanalyzedbyBowtie2,FastQC,MACS,deepTools,phastCons,Rlanguage,HOMERandGeneOntology(GO).Theenrichmentsitesofopenchromatinregionswerehighlyconserved.Althoughtheopenregionsweremainlyenrichedbetweenthepro-moterandtranscriptionstartsite(TSS),theenrichmentsintheyoungandtheoldarenotthesame,varyingwithagingof2BScells.Furthermore,GOanalysiswasperformedtoanalyzethedifferentiallyexpressedgenesindifferentagesof2BScells.Thegenesweremainlyclusteredincellularprocess,metabolicprocess,biologicalregulation,bindingandsoon.Theseresultsshowedthattheagingprocesswasaccompaniedbydeclineddegreeofchromatinopennessanddecreasedleveloftranscriptionalregulation,whichrevealedthecloserelationshipbetweenthereplicativesenescenceandtranscriptionalregulation.Keyword：replicativesenescence;ATAC-seq;transcriptionalregulation;openingofchromatin;復制性衰老又稱細胞衰老,當細胞周期停滯、凋亡抵抗以及部分基因如P16INK4a表示出改變時,細胞出現衰老表型[1,2]。復制性衰老是整體衰老的基礎,最近研究發現其具有拮抗癌癥、使癌細胞得以去除的能力[3,4]。近年來,隨著人口老齡化的發展及癌癥等老年病發病率的增加,細胞衰老作為老年疾病的主要病理基礎,已成為分子生物學研究的熱門之一[5,6,7,8]。在真核生物中,DNA圍繞組蛋白分級折疊,構成核小體,而核小體則通過進一步嚴密折疊構成高級構造染色質[9]。此時,根據細胞所處的狀態,其基因組可劃分為活潑踴躍表示出的基因與相對沉默的基因。在DNA分級折疊的這一經過中,非活潑踴躍基因表示出相關的構造區域相對封閉、核小體構造致密,呈轉錄抑制狀態;而與活潑踴躍表示出基因相關的區域則呈開放狀態,核小體相對松懈,有利于轉錄因子、組蛋白等調節元件與這些基因的啟動子、加強子互相作用,促進轉錄的進行[9,10]。通常情況下,染色質的開放程度可代表轉錄因子、組蛋白等調節元件與DNA的結合狀態,因而根據染色質的開放程度明確生理學狀態、疾病發生發展中表觀遺傳學狀態的改變,也已成為表觀遺傳學技術的主要切入點[9,11,12,13]。迄今為止,細胞衰老經過中表觀遺傳學狀態的改變尚不明確,從染色質空間構造、基因組學角度對染色質開放性的變化進行分析討論的研究較少,在細胞衰老經過中染色質開放程度的詳細變化仍未知[14,15]。故本研究采用基于轉座酶Tn5和高通量測序的染色質開放性分析技術ATAC-seq(assayfortransposase-accessiblechromatinwithhighthroughputsequencing),研究在復制性衰老經過中全基因組染色質開放程度的變化,為探尋求索復制性衰老經過中表觀遺傳以及轉錄調控水平的變化提供初步的研究基礎。1、材料和方式方法1.1、材料年輕代齡(PD26)和年老代齡(PD55)的人胚肺二倍體成纖維細胞(2BS)購自武漢大學中國典型培養物保藏中心。1640培養基、胎牛血清(fetalbovineserum,FBS)、青霉素鏈霉素混合液和磷酸鹽緩沖液(phosphatebuffersaline,PBS)購自美國Gibco公司;Tris、NaCl以及MgCl2購自BBI生命科學有限公司;Tn5轉座酶源自美國Illumina公司;NP-40和nuclease-freeH2O購自上海碧云天生物技術有限公司;DNA純化試劑盒源自美國Qiagen公司;PCR擴增試劑盒源自美國NewEnglandBioLabs公司。二代測序儀(150PE,Illumina公司,美國);qPCR儀(480,Roche公司,瑞士);低溫高速離心機(CR3I,Thermo公司,美國)。1.2、Tn5酶切反響參照文獻[16]的方式方法,對2BS細胞進行Tn5酶切反響,主要步驟如下:取凍存于液氮中的PD26、PD55的2BS細胞進行復蘇。復蘇后的細胞置于T75培養瓶,采用添加10%胎牛血清、1%青霉素鏈霉素混合液的1640培養基培養,培養箱條件為37℃、5%CO2。當細胞長滿且生長狀態良好時,用胰酶將細胞消化,并使用臺盼藍染色,過濾除去死細胞后進行細胞計數,計數50000個細胞。取計數后的細胞,在4℃條件下500g離心5min后棄上清;再參加50L預冷PBS吹打洗滌,4℃條件下500g離心5min,棄上清;參加50L預冷的裂解液(由10mmol/LTrispH7.5、10mmol/LNaCl、3mmol/LMgCl2以及0.1%NP-40配置而成)懸浮細胞,在4℃條件下500g離心10min后棄上清。此時向細胞沉淀中參加50L轉座反響體系(表1),吹打混勻后,37℃孵育30min。待孵育結束后使用DNA純化試劑盒純化DNA,純化完成后使用10L試劑盒中的洗脫液將其洗脫,隨后進行PCR擴增。PCR反響體系如表2所示,循環條件見表3。華而不實,PD26對應的primer1序列為AGGCAGAA,PD55對應的primer1序列為TCCTGAGC;barcodedprimer2序列為AATGAT-ACGGCGACCACCGAGATCTACACTCGTCGGCA-GCGTCAGATGTG。PCR結束后,再次使用DNA純化試劑盒純化DNA[17]。表1轉座反響體系表2PCR反響體系1.3、ATAC-seq測序、質量評估及原始數據過濾使用IlluminaHiseq150PE測序儀對Tn5酶切反響得到的DNA進行雙端測序,辨別chasity0.6的堿基,對每個測序序列(sequencedread)進行passfilter檢驗(當某一read的前25個堿基中,僅存在一個或者不存在chasity0.6的堿基時,斷定該read符合檢驗標準),過檢后的單克隆reads則使用FastQC進行測序質量評估。斷定評估合格的標準如下:每個堿基的質量均5、每組堿基的中位數均tile的質量檢測中不存在暖色結果;序列質量評分均集中于高分;GC堿基含量分布均勻;GC堿基占比均勻;測序中序列長度分布較一致;測序序列中接頭(adapter)出現頻率較低;測序序列中特征性短序列重復出現的次數較低。假如評估合格,則進行數據過濾,去除低質量(QPhred20的堿基數占整個read長度50%以上的reads)、N(N表示無法確定堿基信息)的比例10%和包含adapter的reads,以保證分析質量。表3PCR循環條件1.4、ATAC-seq的序列比對和結合位點(peak)掃描使用Bowtie2(版本2.25)對參考基因組(hg19)建立索引回帖后,進行序列比對。然后根據泊松分布,使用MACS(版本2.1.0)對測序后各個區域基于唯一比對的read數的P-value值進行計算,當P-value1.0E-05時,該區域則被以為是一個結合位點(peak)。標準化后的P-value值越高,理論上該區域的染色質越致密[18]。1.5、保守性分析使用phastCons軟件,根據phylo-HMM模型對結合位點的保守性進行計算[19]。1.6、回帖序列在基因附近的信號分布情況分析將所有基因作為目的區域,使用deepTools對基因長度(轉錄起始位點TSS到轉錄終止位點TTS之間的距離)進行歸一化處理,繪制ATAC-seqreads(使用bigwig文件)在轉錄起始位點上游3kb到轉錄終止位點下游3kb范圍內的平均信號值,以averageplot的形式將結果展示。1.7、結合位點注釋分析使用HOMER軟件的annotatePeaks.pl工具尋找離每個結合位點距離近期的轉錄起始位點所屬的基因(下面稱為結合位點鄰近基因),以及結合位點覆蓋的基因,并注釋其覆蓋的功能性區域。注釋的功能性區域包括啟動子到轉錄起始位點的區域(promoter-TSS,默認范圍為基因起點上游1kb到基因起點下游100bp)、轉錄終止位點(TTS,默認范圍為基因終點上游100bp到基因終點下游1kb)、外顯子(exon)、內含子(intron)和基因間區(intergenic)。假如一些功能性區域有重疊部分,則按上述順序優先選擇排在前面的功能性區域作為最終注釋。對上述注釋分析結果中結合位點在所覆蓋基因的功能性區域上的分布進行統計,并使用R語言的ggplot2包將統計結果繪制成餅狀圖和柱狀圖。1.8、差異結合位點分析使用bedtools工具將2BS-PD55與2BS-PD26的peak文件合并,并對合并文件進行處理。假如兩個結合位點有重疊區域,則合并成一個新的結合位點。計算PD26與PD55在每個合并的新結合位點區域上的read數目,然后使用DEGseq進行差異分析,得到樣本間的差異結合位點,即差異染色質開放區域。本次差異結合位點的挑選原則為|log2FC|1且P-value0.05。1.9、差異結合位點注釋分析使用HOMER軟件的annotatePeaks.pl工具尋找離PD26與PD55差異結合位點距離近期的轉錄起始位點所屬的基因,以及差異結合位點覆蓋的基因,并注釋其覆蓋的功能性區域。對差異結合位點注釋分析結果中結合位點在所覆蓋基因的功能性區域上的分布進行統計,并使用R語言ggplot2包將統計結果繪制成餅狀圖。1.10、GO富集分析基因本體論(GeneOntology,GO)數據庫可對基因功能進行富集分析。我們使用GOTermFinder,將PD26與PD55差異結合位點鄰近的基因與GOterm數據庫比對,計算每個term的基因數量,然后應用超幾何檢驗,找出整個基因組背景下結合位點相關基因中顯著性富集的GOterm,進而辨別出差異結合位點相關基因行使的主要生物學功能。2、結果2.1、測序質量2.1.1、測序的錯誤分布率對Tn5酶切后的年輕與年老細胞中的DNA進行雙端測序,排除酶切后堿基質量因素以及測序儀和測序試劑等的影響后,所得reads的測序錯誤率如此圖1所示,均低于0.06%。本實驗測序的錯誤分布率檢查合格、測序質量較高,提示此次測序結果可用于后續分析。Fig.1Errorratedistributionalongr2.1.2、堿基含量分布使用FastQC對Tn5酶切后的年輕與年老細胞中的DNA進行堿基組成檢測,結果如此圖2所示。各嘌呤堿基和嘧啶堿基的組成比例都較平衡,低質量堿基(指在某一read的某一位置中,含量小于20%的堿基)比例較低且僅出現于隨機引物擴增的開場階段,表示清楚本次測序結果較好,可信度較高。2.2、酶切片段長度的分布酶切片段(insertsize,也稱作fragmentlength)是指經轉座酶Tn5切割后獲得的DNA片段,常用于監測Tn5酶切反響中酶的用量能否適宜。本實驗的酶切片段長度分布如此圖3所示,從左至右依次為無核小體片段(nucleosomefree,即開放染色質區域)、單核小體片段(包括一個核小體的酶切片段,長度一般在147bp到147*2bp之間)和雙核小體片段(包括兩個核小體的酶切片段)的長度分布,表示清楚本次實驗酶切效果較好。2.3、結合位點統計對測序產生的reads的P-value值進行計算,通常標準化后的P-value值越小,理論上該區域的染色質越疏松[18],越可能成為染色質的開放區域進而結合轉錄因子。對本次測序樣本的peaks進行統計分析后發現,PD26的2BS細胞有45681個peaks,而PD55的2BS細胞有50075個peaks。2.4、結合位點保守性分析圖4展示了年輕細胞與年老細胞進化保守性的分析結果。在脊椎動物基因組學中,調節復制、轉錄和翻譯的序列,常被以為是功能性保守序列[20]。在本次實驗中,我們對結合位點進行統計分析發現,年輕細胞與年老細胞的染色質開放性變化區域較保守、類似性高,因而能夠初步推斷這些結合位點的相關區域在轉錄表示出中具有較高的重要性。圖2測序序列堿基分布示意圖Fig.2Sequencecontentacrossallbases圖3酶切片段長度的分布Fig.3FragmentlengthdistributionX軸代表Tn5酶切后DNA片段的長度,Y軸代表該長度片段出現的比例。TheX-axisrepresentsfragmentlength,andtheY-axisrepresentsthefrequencyoflengthdistribution.圖4結合位點保守性分析Fig.4TheresultsofconservationX軸代表此次檢測所得的結合位點對應基因組中該位點的相對位置信息,華而不實0表示相對位置的中點;Y軸代表保守性的評分,分值越高代表保守性越強。TheX-axisrepresentspositioninformation,and0representsrelativepositionofpeaks.TheY-axisrepresentsconservationscore(ahigherscoreindicatesahigherconservation).2.5、回帖序列(mappedreads)在基因附近的平均信號分析在所有基因長度歸一化后,將reads與hg19基因組比對,計算其在基因組上的分布情況,結果如此圖5所示。從結果可知,回帖序列(即染色質的開放區域)主要集中在轉錄起始位點附近,提示這些染色質的開放區域與活潑踴躍的轉錄調控相關。同時,我們發現,隨著細胞衰老程度的加深,轉錄起始位點附近的染色質開放程度明顯下降,提示隨著細胞衰老的發生發展,轉錄水平進一步下降。圖5回帖序列在基因附近的平均信號分布Fig.5TherelationshipbetweenmappedreadsandgenesinaverageplotX軸表示歸一化后的基因范圍,-3.0表示TSS上游3kb,3.0表示TTS下游3kb;Y軸表示位點的平均信號,數值越大代表越富集。TheX-axisrepresentsthenormalizedgenerange,-3.0represents3kbupstreamtheTSS,and3.0represents3kbdownstreamtheTTS.TheY-axisrepresentstheaveragesignalvalueofthesite.TheenrichmentdegreeisproportionaltotheY-axisvalue.2.6、結合位點在全基因組功能區上的注釋全基因組的功能區域分為啟動子到轉錄起始位點的區域、轉錄終止位點、外顯子區、內含子區和基因間區。轉錄因子一般集中分布在轉錄起始位置附近,即TSS區域。結合位點在全基因組功能性區域上的分布如此圖6所示,年輕細胞(2BS-PD26)、年老細胞(2BS-PD26)的染色質開放區段主要位于:內含子區、基因間區、啟動子到轉錄起始位點的區域、外顯子區和終止子。除去在全基因組中占比極大且不介入基因表示出的內含子和基因間區后,染色質開放區主要富集于啟動子到轉錄起始位點的區域。除此之外,圖7結果顯示年輕細胞與年老細胞中結合位點在該區域分布的比例分別為23.04%和15.4%。以上信息表示清楚轉錄調控對復制性衰老有著重要影響,且隨著細胞的逐步衰老,啟動子到轉錄起始位點區域的開放比例進一步下降,即轉錄的水平進一步下降。2.7差異結合位點統計通過分析PD55與PD26之間的差異結合位點,能夠得到衰老狀態下特異性染色質開放區域。本次實驗樣本間的差異結合位點數目為15812。圖6結合位點在基因組功能性區域的分布統計圖Fig.6Theannotationofpeaksinthegenome圖7結合位點在啟動子到轉錄起始位點區域的分布統計圖Fig.7Theannotationofpeaksinpromoter-TSS2.8、差異結合位點在全基因組功能區上的注釋2BS-PD26和2BS-PD55的差異結合位點在全基因組功能區的注釋如此圖8所示。結果表示清楚,細胞衰老后染色質關閉與開放的區域主要集中在內含子區和啟動子到轉錄起始位點區域。華而不實,染色質關閉的區域有40.31%分布于內含子區,有26.45%分布于啟動子到轉錄起始位點區域;染色質開放的區域有52.55%分布于內含子區,有6.18%分布于啟動子到轉錄起始位點區域。同時,在注釋中我們能夠發現基因間區占比擬多,揣測可能由于基因間區在基因組上的占比擬大,故檢測到的結合位點相對其他區域更多,但此類結合位點并非真正的轉錄調控因子的結合位點。圖8差異結合位點在基因組功能性區域的分布統計圖Fig.8Theannotationofdifferentpeaksinthegenome(A)細胞衰老后,染色質關閉區域在基因組功能性區域的分布;(B)細胞衰老后,染色質開放區域在基因組功能性區域的分布。(A)Thedistributionofclosingchromatinregionsannotation;(B)Thedistributionofopenchromatinregionsannotation.2.9、差異位點靶基因的GO分析對2BS-PD26和2BS-PD55的每個差異結合位點的鄰近基因進行GO注釋,部分結果如表4所示。表4差異位點靶基因GO富集分析圖9差異位點靶基因的GO富集分析Fig.9GOtermenrichmentanalysistothetargetgeneofdifferentpeaks選取差異位點鄰近基因顯著富集的GO條目繪制柱狀圖,結果如此圖9所示。富集程度較高的有細胞進程(cellularprocess)、代謝(metabolicprocess)、生物性調節(biologicalregulation)、結合功能(binding)、催化活性(catalyticactivity)等,表示清楚復制性衰老與細胞本身調節嚴密相關,而這些生物學經過,均與轉錄調控有密切聯絡。3、討論盡管當前針對復制性衰老經過中表觀遺傳學變化的研究方式方法有很多,如使用DNA酶1酶切基因組后進行高通量測序的DNase-seq、使用微球菌核酸酶酶切基因組后進行高通量測序的MNase-seq和利用甲醛對染色質固定后進行高通量測序的FAIRE-seq,但此類方式方法對實驗組細胞數量要求極高、所需測序深度較高、實驗步驟繁瑣且可重復性極低[10,16]。ATAC-seq的開發與應用,極大程度減少了樣本的使用量、增加了實驗的可重復性[21]。本實驗中的測序質量分析顯示:堿基錯配率極低、堿基含量分布均勻,這再次證明該實驗方式方法的準確性。同時,ATAC-seq作為檢測表觀遺傳學變化的一種新方式,可對靶基因的表觀狀態給出更為準確的揣測,具有一定的可挖掘性[22,23]。本研究通過ATAC-seq技術結合生物信息學分析方式方法,對年輕代齡(PD26)和年老代齡(PD55)的2BS細胞展開了分析,初步發現,2BS細胞中結合位點的基因保守性均較強;ATAC-seq的回帖序列主要富集于轉錄起始位點附近;結合位點特異性地注釋于啟動子到轉錄起始位點區域,提示復制性衰老與轉錄調控具有較強的相關性。文中在分析回帖序列在基因附近的平均信號分布時,發現染色質的開放區域主要富集在轉錄起始位點附近,且年輕細胞與衰老細胞在轉錄起始位點的開放程度不同,衰老時染色質的開放程度明顯下降。這可能與衰老時基因的表示出水平改變相關,衰老時細胞內各類蛋白質的表示出大部分下降。本實驗室前期研究發現,E2F1、POLD1等基因的表示出水平隨著細胞衰老出現下降的趨勢[24,25,26],此類基因的啟動子區在年輕細胞中呈開放狀態,當細胞衰老時則表現得較為封閉。故此類基因的回帖序列經過富集分析,可出現衰老時轉錄起始位點的富集減少。但是,轉錄起始位點區域富集減少并非意味著衰老時細胞不出現轉錄調控。此時,細胞內以E2F1、POLD1為代表的一類基因還能夠出現微弱表示出,而另外一些基因(如P16INK4a等)則隨著細胞衰老的出現表示出上升[27,28],這些表示出上調基因的啟動子區會隨衰老的發生發展進一步開放,進而使轉錄起始位點富集的回帖序列增加。由于衰老時呈表示出增加的基因仍較表示出減少的基因少,所以從全基因組的層面分析,衰老時染色質的開放區域在轉錄起始位點區的富集下降。其實,結合位點在全基因組功能區上的注釋,也進一步說明了這一問題。對于年輕與年老的2BS細胞,其開放區段均富集于啟動子到轉錄起始位點區域,表示清楚轉錄調控對復制性衰老有著重要影響。隨著細胞的逐步衰老,啟動子區的開放比例進一步下降,講明當細胞衰老發生時,細胞內的轉錄調控水平減弱,而在年輕細胞中轉錄調控較為活潑踴躍。當然,在年老細胞中也存在開放的啟動子區,這是由于衰老的細胞也存在需要轉錄表示出的基因。除此之外,我們使用GO富集分析進一步驗證了轉錄調控對復制性衰老的調節。華而不實,牽涉細胞進程、代謝、生物性調節、結合功能、催化活性等功能的基因,在細胞衰老經過中富集,此類基因與細胞本身調節、轉錄調控均有密切聯絡。同時,我們也發現,作為復制性衰老經典標志物的P53、P16INK4a[1,29,30]等,均在GO分析結果中富集,如P53富集于GO:0048522、GO:0071840、GO:0016043等;P16INK4a富集于GO:0071840、GO:0016043、GO:0044260等。綜上所述,從轉錄調控探究、理解復制性衰老,應為我們后續研究的主要工作之一。以下為參考文獻[1]LOPEZ-OTINC,BLASCOMA,PARTRIDGEL,etal.Thehallmarksofaging[J].Cell,2020,153(6):1194-1217.[2]CHILDSBG,BAKERDJ,KIRKLANDJL,etal.Senescenceandapoptosis:duelingorcomplementarycellfates?[J].EMBOReports,2020,15(11):1139-1153.[3]CAMPISIJ,DADDADFF.Cellularsenescence:whenbadthingshappentogoodcells[J].NatureReviewsMolecularCellBiology,2007,8(9):729-740.[4]CHILDSBG,DURIKM,BAKERDJ,etal.Cellularsenescenceinagingandage-relateddisease:frommechanismstotherapy[J].NatureMedicine,2021,21(12):1424-1435.[5]PALS,TYLERJK.Epigeneticsandaging[J].ScienceAdvances,2021,2(7):e1600584.[6]KENNEDYBK,BERGERSL,BRUNETA,etal.Geroscience:linkingagingtochronicdisease[J].Cell,2020,159(4):709-713.[7]MOSKALEVAA,ALIPERAM,SMIT-MCBRIDEZ,etal.Geneticsandepigeneticsofagingandlongevity[J].CellCycle,2020,13(7):1063-1077.[8]CHANDRASEKARANA,IDELCHIKMDPS,MELENDEZJA.Redoxcontrolofsenescenceandage-relateddisease[J].RedoxBiology,2021,11:91-102.[9]VOSSTC,HAGERGL.Dynamicregulationoftranscriptionalstatesbychromatinandtranscriptionfactors[J].NatureReviewsGenetics,2020,15(2):69-81.[10]BUENROSTROJD,GIRESIPG,ZABALC,etal.Transpositionofnativechromatinforfastandsensitiveepigenomicprofilingofopenchromatin,DNA-bindingproteinsandnucleosomeposition[J].NatureMethods,2020,10(12):1213-1218.[11]LIB,CAREYM,WORKMANJL.Theroleofchromatinduringtranscription[J].Cell,2007,128(4):707-719.[12]RUIZJL,TENAJJ,BANCELLSC,etal.CharacterizationoftheaccessiblegenomeinthehumanmalariaparasitePlasmodiumfalciparum[J].NucleicAcidsResearch,2021,46(18):9414-9431.[13]TSOMPANAM,BUCKMJ.Chromatinaccessibility:awindowintothegenome[J].EpigeneticsChromatin,2020,7:33.[14]QUANH,YANGY,LIUS,etal.ChromatinstructurechangesduringvariousprocessesfromaDNAsequenceview[J].CurrentOpinioninStructuralBiology,2022,62:1-8.[15]YUR,MCCAULEYB,DANGW.Lossofchromatinstructuralintegrityisasourceofstressduringaging[J].HumanGenetics,2020,139(3):371-380.[16]BUENROSTROJD,WUB,CHANGHY,etal.ATAC-seq:amethodforassayingchromatinaccessibilitygenome-wide[J].CurrentProtocolsinMolecularBiology,2021,109:21-29.[17]CORCESMR,TREVINOAE,HAMILTONEG,etal.AnimprovedATAC-seqprotocolreducesbackgroundandenablesinterrogationoffrozentissues[J].NatureMethods,2021,14(10):959-962.[18]ZHENGR,WANC,MEIS,etal.CistromeDataBrowser:expandeddatasetsandnewtoolsforgeneregulatoryanalysis[J].NucleicAcidsResearch,2022,47(D1):D729-D735.[19]FANX,ZHUJ,SCHADTEE,etal.StatisticalpowerofphyloHMMforevolutionarilyconservedelementdetection[J].BioMedCentralBioinformatics,2007,8:374.[20]SIEPELA,BEJERANOG,PEDERSENJS,etal.Evolutionarilyconservedelementsinvertebrate,insect,worm,andyeastgenomes[J].GenomeResearch,2005,15(8