2017-01-22 43 views
3

因此,從HumanGenome19項目(hg19)的矩陣中,我想要生成一個圖,其中我將繪製基因和它們的鏈。在R中爲HG19生成基因圖

hg19的數據看起來像這樣:

GENE CHR txStart txEnd Size STRAND 
RBBP8 chr18 2050000 42016610 113940 - 
CCND3 chr18 41902670 42016610 113940 - 
GGNBP1 chr18 33551475 33556803 5328 + 
LINC00336 chr18 33553882 33561115 7233 - 
PGM3 chr18 83874592 83903655 29063 - 
PGM3 chr18 83874592 83903012 28420 - 
PGM3 chr18 83874592 83903012 28420 - 

+指該基因將在right方向被繪製,而 - 表示left

你知道任何R軟件包可以產生如下給出的HG數據嗎?

另外,我不關心基因在Y-axis上的位置,它們應該定位,這樣名字纔會清晰可讀。

enter image description here

免責聲明:劇情,並沒有我所列出的數據作出反應。

+1

箭頭的長度是否意味着對應的東西?像基因大小? – JustGettinStarted

+1

看看http://www.tengfei.name/ggbio/ ...也許Capter 7的小插曲可能會給你一個提示 – Drey

+0

@JustGettinStarted是的。 – Jack

回答

0

我無法弄清楚如何爲您的數據做到這一點,因爲範圍非常寬。我使用arrows製作了類似的假數據圖。但是,如果沒有進一步修改,這可能無法正常工作。

DATA

mydata = structure(list(GENE = structure(c(5L, 1L, 2L, 3L, 4L, 4L, 4L), .Label = c("CCND3", 
"GGNBP1", "LINC00336", "PGM3", "RBBP8"), class = "factor"), CHR = structure(c(1L, 
1L, 1L, 1L, 1L, 1L, 1L), .Label = "chr18", class = "factor"), 
    txStart = c(20500000L, 20780190L, 20982780L, 21218290L, 21533530L, 
    21851180L, 22073300L), txEnd = c(20557770L, 20806930L, 21140420L, 
    21299010L, 22513330L, 21863505L, 22162610L), Size = c(57770L, 
    26740L, 157640L, 80720L, 979800L, 12325L, 89310L), STRAND = structure(c(1L, 
    1L, 2L, 1L, 1L, 1L, 1L), .Label = c("-", "+"), class = "factor")), .Names = c("GENE", 
"CHR", "txStart", "txEnd", "Size", "STRAND"), class = "data.frame", row.names = c(NA, 
-7L)) 

CODE

graphics.off() 
windows(width = 10, height = 7) 
plot(x = c(min(mydata$txStart), max(mydata$txEnd)), y = c(0,nrow(mydata)+1), 
    type = "p", pch = NA, axes = FALSE, xlab = "Chromosome ## Mb", ylab = "") 

# Divide range of genes into 4 groups to add ticks and labels in next step 
x_label = seq(min(mydata$txStart),max(mydata$txEnd),(max(mydata$txEnd) - min(mydata$txStart))/4) 

#add labels at 4 different places on x-axis 
axis(1, at = x_label, labels = paste(as.character(round(x_label/1000000),2)," Mb",sep ="")) 

y_pos = 1 #Starting vertical position of genes 
for (i in 1:nrow(mydata)){ 
    #use arrows from txStart to txEnd. Based on STRAND value, this may have to change with if else 
    arrows(x0 = mydata$txStart[i], x1 = mydata$txEnd[i], y0 = y_pos, y1 = y_pos, length = 0.05) 

    #Obtain x position to put gene label 
    x_pos = mydata$txStart[i] + (mydata$txEnd[i] - mydata$txStart[i])/2 
    gene_label = paste(mydata$GENE[i]) 

    #Add gene label 
    text(x_pos, y_pos, bquote(italic(.(gene_label))), pos = 3, col = "darkgrey", cex = 0.8) 
    y_pos = y_pos + 1 #Comment this if you want all genes on the same level 
} 

PLOT enter image description here