최근에 종양 내과에 계신 선배와 함께 담관암 (Biliary tract cancer) 환자들의 암유전체 (Cancer Genomics) 관련 NGS 분석을 시작했습니다. 확실히 작년부터 다양한 질환과 환자들의 다양한 유전체 데이터를 접하다보니, 데이터셋의 특성에 따라서 분석 및 접근 방법이 많이 다른 것을 느낍니다. 사실 그동안 저는 주로 Germline variant 분석을 했었는데, Somatic variant 가 더 중요한 암 환자들을 분석하기 위해서는 추가로 더 공부하고 알아야 할 내용들이 많은 것 같습니다. 특히 암과 같은 경우는 선천적으로 가지고 있는 Germline variant와 살면서 축적된 Somatic mutation을 종합적으로 함께 고려해야하기 때문에, 더 복잡한 측면이 있는 것 같습니다.
관련 포스팅 보기 > [유전학 중요개념 정리] Germline vs. Somatic mutation
오늘은 암유전체 분석 관련 논문들의 Figure 1을 차지하는 Waterfall plot을 만드는 법을 잘 소개하고 있는 페이지가 있어 관련 내용을 스크랩 합니다.
Waterfall plot: 환자군에서 나타나는 Mutation을 유전자별로 나타내어, 전체 암 유전체의 특성을 잘 나타내주는 plot으로 마치 폭포가 떨어지는 모양과 비슷한데서 그 이름의 유래가 있습니다.
관련 설명 링크 가기:
Introduction to waterfall plots (Griffith Lab)
이제 연구자들에게 Bioconductor와 R과 같은 프로그램은 필수죠. 다행히 저와 같은 코드를 잘 모르는 사람을 위해서 R package가 잘 만들어져 있습니다. 위 페이지에서 소개하는 필수 코드를 약간 정리하여 올립니다.
library(GenVisR) setwd("c:/BTC_R/sample") # Load data file mutationData <- read.csv("BKM120_mutationdata.csv") clinicalData <- read.csv("BKM120_clinical.csv") mutationBurden <- read.csv("BKM120_mutationburden.csv") # Reformat the mutation data for waterfall() mutationData <- mutationData[,c("patient", "gene.name", "trv.type", "amino.acid.change")] colnames(mutationData) <- c("sample", "gene", "variant_class", "amino.acid.change") # Create a vector to save mutation priority order for plotting mutation_priority <- as.character(unique(mutationData$variant_class)) # Create an initial plot waterfall(mutationData, fileType = "Custom", variant_class_order=mutation_priority) # Define a mutation hierarchy mutationHierarchy <- c("nonsense", "frame_shift_del", "frame_shift_ins", "in_frame_del", "splice_site_del", "splice_site", "missense", "splice_region", "rna") # define colours for all mutations mutationColours <- c("nonsense"='#4f00A8', "frame_shift_del"='#A80100', "frame_shift_ins"='#CF5A59', "in_frame_del"='#ff9b34', "splice_site_del"='#750054', "splice_site"='#A80079', "missense"='#009933', "splice_region"='#ca66ae', "rna"='#888811') # Find which samples are not in the mutationBurden data frame # First, let's look at the sample names in the mutationData and mutationBurden mutationData$sample mutationBurden$sample # Now, determine the non-overlap between these values sampleVec <- unique(mutationData$sample) sampleVec[!sampleVec %in% mutationBurden$sample] # Fix mutationBurden to match mutationData mutationBurden$sample <- gsub("^WU(0)+", "", mutationBurden$sample) # Check for non-overlap again sampleVec[!sampleVec %in% mutationBurden$sample] # reformat clinical data to long format library(reshape2) clinicalData_2 <- clinicalData[,c(1,2,3,5)] colnames(clinicalData_2) <- c("sample", "Months on Study", "Best Response", "Treatment Setting") clinicalData_2 <- melt(data=clinicalData_2, id.vars=c("sample")) # find which samples are not in the mutationBurden data frame sampleVec <- unique(mutationData$sample) sampleVec[!sampleVec %in% clinicalData$sample] # fix mutationBurden to match mutationData clinicalData_2$sample <- gsub("^WU(0)+", "", clinicalData_2$sample) # create the waterfall plot waterfall(mutationData, fileType = "Custom", variant_class_order=mutationHierarchy, mainPalette=mutationColours, mutBurden=mutationBurden, clinData=clinicalData_2, clinLegCol=3, clinVarCol=c('0-6'='#ccbadc', '6.1-12'='#9975b9', '12.1+'='#663096', 'Partial Response'='#c2ed67', 'Progressive Disease'='#E63A27', 'Stable Disease'='#e69127', '1'='#90ddee', '2'='#649aa6', '3+'='#486e77'), clinVarOrder=c('1', '2', '3+', 'Partial Response', 'Stable Disease', 'Progressive Disease', '0-6', '6.1-12', '12.1+'), section_heights=c(1, 5, 1), mainLabelCol="amino.acid.change", mainLabelSize = 3) # Create a sample ordering sample_ordering <- c("19", "5", "31", "22", "12", "2", "32", "8", "28", "18", "4", "24", "23", "17", "11", "14") # Create a gene ordering gene_ordering <- c("CDH1", "MALAT1", "RUNX1", "NCOR1", "GATA3", "FOXA1", "ESR1", "CBFB", "TBX3", "TAB1", "MED12", "XBP1", "TP53", "RB1CC1", "BRCA2", "ATM", "SMARCD1", "MLL3", "MLL2", "ARID1A", "FBXW7", "CAV1", "MAP3K1", "MAP2K4", "NOTCH4", "PDGFRA", "ERBB3", "ERBB2", "RELN", "MAGI3", "MTOR", "AKT2", "AKT1", "PTEN", "PIK3CA") # Create a gene ordering waterfall(mutationData, fileType = "Custom", variant_class_order=mutationHierarchy, mainPalette=mutationColours, mutBurden=mutationBurden, clinData=clinicalData_2, clinLegCol=3, clinVarCol=c('0-6'='#ccbadc', '6.1-12'='#9975b9', '12.1+'='#663096', 'Partial Response'='#c2ed67', 'Progressive Disease'='#E63A27', 'Stable Disease'='#e69127', '1'='#90ddee', '2'='#649aa6', '3+'='#486e77'), clinVarOrder=c('1', '2', '3+', 'Partial Response', 'Stable Disease', 'Progressive Disease', '0-6', '6.1-12', '12.1+'), section_heights=c(1, 5, 1), mainLabelCol="amino.acid.change", mainLabelSize=3, sampOrder=sample_ordering, geneOrder=gene_ordering)