P值與基因組學（1）：從fastq文件的分析的分析談起

計劃寫一個系列，可能到十月中旬就寫完了。

測序儀在每測一個鹼基的時候會進行一次假設檢驗，算出如果該鹼基是由於噪音而測出的概率p，更一般的情況下（在後面的文章我們會提到的），我們稱該值為p值

通常情況下，p是很小的，對付很小的書和很大的數，我們有對數這個神器，於是引入了Q socre:

(1)

我們可以發現，p=0.1,Q=10; p=0.01,Q=20,p=0.001,Q=30，即：如果該鹼基是由於噪音測出的概率越小，p越小，Q越大，所以Q值越大，該鹼基的質量越好。

比如圖中的鹼基T的Q score 為25， 根據換公式(1)，我們得到在該位置由於雜訊測得鹼基T的概率為。

現在問題來了，我們如何衡量如圖中所示的這樣一段讀長的序列的質量呢？

為了深入討論，我們先講一個思路，

假設我們對現在所得的一條Reads,序列讀長為,位置1的鹼基的p值為,位置2的鹼基的p值為，，，位置的鹼基的p值為,

假設其中最大的一個p值為，對於任意的,我們有

那麼對於任意位置，該位置的準確度為，我們有：

記該序列的p值為，假設之間是相互獨立的，則我們有：

(2)

根據（2），我們如果取n=72，,則我們有

誒，巧了，小於費歇爾老爺子說的0.05。

但是問題來了，隨著測序技術的不斷發展，這個讀長是在不斷增長的，所以根據（2）當n越大的時候，遲早會變得大於0.05，更糟糕的是，對n取極限的話，的上界會變為1，這真是一個悲傷的故事，那麼我們該怎麼辦。

我們可以用 Bonferroni correction（Bonferroni_correction）來解決我們的問題，將（2）改寫為：

(3)

如果繼續取,則當我們有n=72時，，n=1000時，由（2）得，，有（3）得到,對n取極限，我們會有

當p甚小時

但是這樣也有個問題，我們如何確定檢驗的閾值呢？如果是Bonferroni correction，比如n=36時，閾值不再是，而是但是我們算出來的p值並不小於0.00139這又讓人不舒服的事情。

而且，對於同一批次的話，我們的Reads不只有一條，而是有很多條的，我們如何確定這一批次的質量呢，且聽下回分解。