序列比對實驗報告
一.實驗內(nèi)容
1.利用序列比對線上工具做序列比較
2.Blast線上工具使用
3.使用r語言實現(xiàn)Needleman-Wunsch算法
二.實驗?zāi)康?/h1>
1.掌握序列比對線上工具使用
2.掌握雙序列比對算法——Needleman-wunsch算法
三.實驗數(shù)據(jù)工具及步驟
1. 利用序列比對線上工具做序列比較
在Swiss-port下載蛋白序列,以BRCC3_HUMAN和 BRCC3_MOUSE的蛋白質(zhì)序列為例编兄,利用EMBL 網(wǎng)站的雙序列比對工具
2.Blast線上工具使用
同樣利用NCBI網(wǎng)站blast在線分析藏姐,blastp
3.實現(xiàn)Needleman-Wunsch算法
1)先把已知的替換積分矩陣導(dǎo)入
替換記分矩陣
2)把要比對的序列文件導(dǎo)入,直接寫兩行序列
例如:
序列
3)實現(xiàn)用Needleman-Wunsch 算法得出打分矩陣蝙茶,根據(jù)公式和替換記分矩陣算出
公式
公式
四.實驗代碼
setwd("F:\\實驗\\轉(zhuǎn)錄組學(xué)\\實驗一")
matrix<-read.table("matrix.txt",header=T) #導(dǎo)入打分矩陣宵呛,行列名AGCT
str(matrix)
colnames(matrix)<-c("A","G","C","T") #設(shè)置列名
替換記分矩陣
seqdata<-read.table(“seqdata.txt”,as.is=T)
seqdata#導(dǎo)入序列茄厘,如圖
序列
#轉(zhuǎn)化成單個字符
seqdata<-as.matrix(seqdata)
seqdata1<-seqdata[1,]#提取序列1
seqdata2<-seqdata[2,]#提取序列2
#統(tǒng)計序列長度
M<-nchar(seqdata1)
N<-nchar(seqdata2)
seqdata1<-strsplit(seqdata1,"",fixed=T)
seqdata2<-strsplit(seqdata2,"",fixed=T)
zseqdata1<-as.character(unlist(seqdata1))
zseqdata2<-as.character(unlist(seqdata2)) #zseqdata1和zseqdata2是轉(zhuǎn)化成單個字符后的序列
#Needleman-Wunsch 算法
gap=-5#已知gap
scorematrix<-matrix(0,N+1,M+1)#構(gòu)造空矩陣,N+1行,M+1列
rownames(scorematrix)<-c(0,zseqdata2)
colnames(scorematrix)<-c(0,zseqdata1)
#計算第一行第一列
scorematrix[1,1]=0
for (i in 0:N+1)
scorematrix[i,1]=gap*(i-1)
for (j in 0:M+1)
scorematrix[1,j]=gap*(j-1)
#計算剩下的
for (i in 1:N+1)
for (j in 1:M+1)
{
scorematrix[i,j]=max(c(scorematrix[i-1,j-1]+matrix[rownames(scorematrix)[i],colnames(scorematrix)[j]],
scorematrix[i-1,j]+gap,
scorematrix[i,j-1]+gap))
}
scorematrix
五.實驗結(jié)果:
分析:Gap open越大,比對空位減少,得分越低,gap越集中
Gap extend變化,比對結(jié)果沒有發(fā)生變化,而gap越分散
PAM-n矩陣稠茂,n越大,序列相似度越低,BLOSUM-n矩陣柠偶,n越大,序列相似度越高
回溯表示R語言代碼現(xiàn)在還沒想出來睬关,如果有寫出來的小伙伴可以交流分享吖
補充:Needleman-wunsch算法原理是設(shè)置打分矩陣诱担,根據(jù)適當(dāng)?shù)拇蚍止絹韺?yīng)的堿基進行打分,有四種情況:1.兩堿基完全匹配2.不匹配3.第一條序列引入空位4.第二條序列引入空位
具體算法:
替換記分矩陣
公式
已知gap=-5
1.寫出替換打分矩陣
依次算出该肴,最終得到替換打分矩陣
2.寫出比對序列
比對結(jié)果:
最終得分為右下角的數(shù)字
score=21
從這開始,依次往回找箭頭藐不,如圖藍色箭頭
書寫比對結(jié)果:先把第一個序列寫出來
A C G T C
然后從最左邊開始寫匀哄,橫箭頭和豎箭頭表示字母對空,斜箭頭表示字母對字母雏蛮,第一個是A對A,第二個是C對空涎嚼,依次對應(yīng),結(jié)果如圖:
image.png
比對結(jié)果.png
斜箭頭代表第一個對應(yīng)第二個
橫箭頭代表第一個對空
豎箭頭代表空對第二個
這樣看來,序列比對四不四挺簡單呢