SFS解釋原文：joint SFS – The G-cat
參考博客：種群動態歷史評估方法之 SFS（Site/Allele Frequency Spectrum） - 簡書 (jianshu.com)
基因組數據統計可以總結數據集中等位基因的變異或分布，而無需我們所有樣品的整個基因序列。我們可以選擇一個非常有效的匯總統計量SFS（等位基因頻譜），這與等位基因頻率測量（如Fst）混淆，SFS本質上是某些等位基因在我們的數據集中的頻率的直方圖。SFS根據每個等位基因的常見程度 將其分類為待定類別，並計算以該頻率出現的等位基因的數量。類別的總數將是可能等位基因的最大數量：對於二倍體物種，每條染色體具有兩個拷貝，這意味着包含的樣本數量是其兩倍。
SFS只是我們可以用來調查種群歷史的眾多工具之一。通過結合基因組技術，聚結理論和更強大的分析方法，SFS似乎准備解決地球上生命進化史的更微妙和更復雜的問題。
- SFS的應用
我們可以對種群歷史進行推斷，而無需模擬大型和復雜的基因序列，而是使用SFS。將我們觀察到的SFS與瓶頸的模擬場景進行比較，並比較預期的SFS，使我們能夠估計該場景的可能性。
例如，我們可以預測，在種群中最近出現遺傳瓶頸的情況下，種群中罕見的等位基因將由於遺傳漂移而不成比例地丟失。因此，SFS的整體形狀將急劇向右移動，留下瓶頸的明確遺傳信號。這與瓶頸的合並測試具有相同的理論背景。
瓶頸如何導致 SFS 偏移的一個代表性示例。中間圖：隨時間變化的等位基因圖，在瓶頸期間丟失了較稀有的變異（黃色和深藍色），和但更常見的變異保留下來（紅色）。左圖：這種趨勢反映在這些等位基因的合並樹上，紅色十字表示該等位基因完全喪失。右圖：所描繪等位基因的瓶頸事件之前（紅色）和之后（藍色）的SFS。在瓶頸之前，變異以通常的指數形式分布：然而，之后，稀有變體的不成比例的損失會導致分布趨於平緩。通常，SFS將由比這里顯示的更多的等位基因構建而成，並且延伸得更遠。
相比之下，一個大型或不斷種群將具有更多來自突然生長和遺傳變異增加的稀有等位基因。因此，與瓶頸相反，SFS分布將偏向頻譜的左端，並具有過多的低頻變異。
與上圖類似，但這次是擴展事件而不是瓶頸。種群的擴張以及隨后的Ne的增加，促進了遺傳漂移中新等位基因的突變（或減少漂移中等位基因的損失），導致更多新的稀有等位基因出現。左圖合並樹，右圖SFS的偏移。
SFS甚至可用於檢測自然選擇下的等位基因。對於基因組的強選擇部分，等位基因應以高（如果選擇為正）或低（如果選擇為負）的頻率發生，並具有更多中間頻率的。
- 評估SFS (estimate Site Allele Frequency)
使用ABGSD軟件進行SFS Estimation，missing data和低depth位點。以bam文件為input。
軟件下載目錄：GitHub - ANGSD/angsd: Program for analysing NGS data.
軟件教程：SFS Estimation - angsd (popgen.dk)
#first generate .saf file
./angsd -bam bam.filelist -doSaf 1 -out smallFolded -anc chimpHg19.fa -GL 2
#now try the EM optimization with 4 threads
misc/realSFS smallFolded.saf.idx -maxIter 100 -P 4 >smallFolded.sfs
#in R
sfs<-scan("smallFolded.sfs")
barplot(sfs[-1])
針對vcf中miss不多情況，可以直接用easySFS軟件
軟件下載：https://github.com/isaacovercast/easySFS
python3 /home/software/easySFS/easySFS.py -i vcf -p ./pops_file.txt -a -f --proj 20,16
第一個群體10個體，第二個群體8個體
最后得到的結果格式.obs
3.Fastsimcoal2，SFS建模模擬種群動態歷史
fastsimcoal安裝和使用教程：Demographic modeling with fastsimcoal2 | Speciation & Population Genomics: a how-to-guide (speciationgenomics.github.io)
使用Fastsimcoal2.6，利用SFS模擬種群動態歷史需要三個輸入文件：
- .tpl 模版文件（參數文件）：用來指定歷史事件、migration、樣本大小、突變率重組率等等。注意，這里的突變率重組率都為【per generation】
- .est 參數評估文件（prior文件）：在這里對參數施加prior，指定參數之間關系等等。
- .obs SFS觀測文件：就是上一步生成的。
這三個文件的前綴必須相同，在這里為"6.GL1.folded.HN"。如果為unfolded SFS，SFS觀測文件后綴應該是DAF之類的而不是MAF，具體得看說明書。
一般來說，有兩種構建模型的方法：
1）可以構建很多模型（可能3～10個）來對種群歷史建模，然后比較每個模型的lnL（也就是LRT檢驗），看那個模型最合適。比如第一個模型是有效群體大小先上升后下降在上升，第二個是先下降在上升等等。對於單個群體（不考慮migration）的簡單模型，一般也就是設置幾個bottleneck時間點，擴張時間點以及對應的種群大小參數。
2）只構建一個模型，限制“時間點”prior的范圍，但是放松Ne prior的范圍，並保持一致，讓算法自己探索種群在某個時間點是擴張還是收縮。在RULES中只對時間做限制，而不對Ne做限制。這種就需要增加運行的次數，讓模型達到收斂，比如設置-L 100。例如：
.est文件：
.tpl文件：
設置好配置文件以后，開跑：
# 對於folded SFS，記得設置 -m
fsc26 -t ${i}.GL1.folded.HN.tpl -e ${i}.GL1.folded.HN.est -m -0 -C 10 -n 100000 -L 40 -s 0 -M -q -c 5 &
網友評價：生成文件里有likelihood的評估，拿過來比較模型就行了。對最優模型可以再跑bootstrap。
我個人覺得對於單個種群而言，用SFS不一定是好的選擇。PSMC從原理上更加靠譜一些。而fastsimcoal模型的優勢在於其可以考慮多個種群的migration、融合與分裂、群體大小變化等等。
Stairway plot v2
分析方法：https://raw.githubusercontent.com/xiaoming-liu/stairway-plot-v2/master/READMEv2.1.pdf
后繼如果看到SFS相關文獻和分析方法，會續上！