BLOG
Amazon EMR에서 Apache Spark와 함께 ADAM 및 Mango를 사용하는 게놈 데이터셋의 추출 데이터 분석하기
게놈 염기 서열 분석 비용이 급격히 줄어들면서 지난 몇년간 공개적으로 이용 가능한 게놈 데이터의 양이 급증했습니다. 새로운 코호트(cohort)와 연구는 100,000명 이상으로 구성된 방대한 데이터셋을 생성했습니다. 동시에, 이러한 데이터셋은 모집단 전체의 유전적 변이를 추출하기 위해 처리되었고, 각 코호트(cohort)에 대한 대량의 변동 데이터를 생성했습니다.
빅 데이터 시대에 Apache Spark와 같은 툴은 대용량 데이터셋의 배치 처리를 위한 사용자 친화적인 플랫폼을 제공했습니다. 하지만 현재의 생물 정보 파이프라인을 충분히 대체할 수 있는 툴을 사용하기 위해서는 유전체 데이터를 처리하기 위해 좀 더 접근하기 쉽고 포괄적인 API가 필요합니다. 또한 이러한 처리된 데이터셋을 양방향으로 탐색하기 위한 지원이 필요합니다.
ADAM 및 Mango는 Apache Spark에서 대규모 게놈 데이터셋을 처리, 필터링 및 시각화할 수 있는 통합 환경을 제공합니다. ADAM을 사용하면 Apache Spark에서 데이터를 집계하고 선택하기 위한 SQL 인터페이스인 Spark SQL을 사용하여 원시 게놈 및 변동 데이터를 프로그래밍 방식으로 로드, 처리 및 선택할 수 있습니다. Mango는 Jupyter 노트북 환경에서 원시 및 집계된 게놈 데이터를 시각화할 수 있도록 지원하므로 큰 데이터셋에서 복수의 해상도로 결론을 도출할 수 있습니다.
ADAM 및 Mango의 결합된 기능을 사용하면 통합환경에서 데이터셋을 로드, 쿼리 및 탐색할 수 있습니다. 단일 노드 생물 정보학 툴을 사용하여 이전에 불가능했던 규모의 게놈 데이터를 쌍방향으로 탐색할 수 있습니다. 이 글에서는 Amazon EMR에서 ADAM 및 Mango를 설정하고 실행하는 방법에 대해 설명합니다. 이러한 툴을 쌍방향적 노트북 환경에서 사용하여 1000 유전체 데이터셋을 탐색하는 방법과 Amazon S3에서 공용 데이터셋으로 공개적으로 이용 가능합니다.
Amazon EMR에서 ADAM 및 Mango 구성
먼저, EMR 클러스터를 시작하고 구성합니다. Mango는 Docker 컨테이너를 사용하여 Amazon EMR에서 쉽게 실행됩니다. 스크립트는 /home/hadoop/mango-scripts에서 사용할 수 있습니다.
aws emr create-cluster
–release-label emr-5.14.0 \
–name ’emr-5.14.0 Mango example’ \
–applications Name=Hadoop Name=Hive Name=Spark \
–ec2-attributes KeyName=<your-ec2-key>,InstanceProfile=EMR_EC2_DefaultRole \
–service-role EMR_DefaultRole \
–instance-groups \ InstanceGroupType=MASTER,InstanceCount=1,InstanceType=c5.4xlarge \ InstanceGroupType=CORE,InstanceCount=4,InstanceType=c5.4xlarge \ –region <your-aws-region> \
–log-uri s3://<your-s3-bucket>/emr-logs/ \
–bootstrap-actions \
Name=’Install Mango’, Path=”s3://aws-bigdata-blog/artifacts/mango-emr/install-bdg-mango-emr5.sh”
Mango 노트북을 시작하려면 다음을 실행합니다.
/home/hadoop/mango-scripts/run-notebook.sh
이 파일은 Amazon EMR의 Docker에서 Mango를 실행하는 데 필요한 모든 환경 변수를 설정합니다. 터미널에서 Mango 노트북 세션의 포트 및 Jupyter 노트북 토큰을 볼 수 있습니다. EMR 클러스터에 대한 마스터 노드의 공용 DNS URL에서 이 포트로 이동합니다.
1000 유전체 프로젝트에서 데이터 로드하기
이제 작업 환경이 구축되었으므로, ADAM 및 Mango를 사용하여 3개의 유전자 배열 데이터(어머니, 아버지 및 아이의 데이터)에서 아이의 흥미로운 변형을 검색할 수 있습니다. 이 데이터는 1000 유전체 프로젝트 AWS 공용 데이터셋에서 사용할 수 있습니다. 이 분석에서는 3개(NA19685, NA19661 및 NA19660)를 보고 부모에게는 없지만 아이에게만 존재하는 변형이 있는지 검색합니다.
특히, 새로운 변형(de novo variants)으로 알려져 있는 부모가 아닌 아이에게서 발견되는 유전자 변형을 확인하고자 합니다. 이것들은 여러가지 질환에 기여할 수 있는 새로운 변형의 광경을 보여 줄 수 있기 때문에 흥미로운 지역입니다.
이러한 예제가 포함된 Jupyter 노트북은 Mango의 GitHub 저장소나 Mango의 실행중인 Docker 컨테이너에 있는 /opt/cgl-docker-lib/mango/example-files/notebooks/aws-1000genomes.ipynb에서 찾을 수 있습니다.
먼저 ADAM 및 Mango 모듈과 필요한 Spark 모듈을 가져옵니다.
# Import ADAM modules
from bdgenomics.adam.adamContext import ADAMContext
from bdgenomics.adam.rdd import AlignmentRecordRDD, CoverageRDD
from bdgenomics.adam.stringency import LENIENT, _toJava
# Import Mango modules
from bdgenomics.mango.rdd import GenomicVizRDD
from bdgenomics.mango.QC import CoverageDistribution
# Import Spark modules
from pyspark.sql import functions as sf
그런 다음 Spark 세션을 만듭니다. 이 세션을 사용하여 버전에 따라 SQL 쿼리를 실행합니다.
# Create ADAM Context
ac = ADAMContext(spark)
genomicRDD = GenomicVizRDD(spark)
Spark SQL을 사용한 변형 분석
염색체 17로부터 변형 데이터의 하위 집합에 로드:
genotypesPath = ‘s3://1000genomes/phase1/analysis_results/integrated_call_sets /ALL.chr17.integrated_phase1_v3.20101123.snps_indels_svs.genotypes.vcf.gz’
genotypes = ac.loadGenotypes(genotypesPath)
# repartition genotypes to balance the load across memory
genotypes_df = genotypes.toDF()
데이터프레임의 열을 인쇄하여 스키마를 확인할 수 있습니다.
# cache genotypes and show the schema
genotypes_df.columns
이 유전자형 데이터셋에는 1000 유전체 프로젝트의 모든 샘플이 포함되어 있습니다. 따라서 다음 번에 유전자형을 필터링하여 NA 19685 트리오에 속하는 샘플만 고려하고 그 결과를 메모리에 캐시합니다.
# trio IDs
IDs = [‘NA19685’, ‘NA19661′,’NA19660’]
# Filter by individuals in the trio
trio_df = genotypes_df.filter(genotypes_df[“sampleId”].isin(IDs))
trio_df.cache()
trio_df.count()
그 다음, 각 변형의 게놈 위치를 결정하는 새로운 열을 데이터프레임에 추가합니다. 이는 염색체(contigName)와 변형의 시작 및 끝 위치로 정의됩니다.
# Add ReferenceRegion column and group by referenceRegion
trios_with_referenceRegion = trio_df.withColumn(‘ReferenceRegion’,
sf.concat(sf.col(‘contigName’),sf.lit(‘:’), sf.col(‘start’), sf.lit(‘-‘), sf.col(‘end’)))
이제 새로운 변형을 찾기 위해 데이터셋을 쿼리할 수 있습니다. 하지만 먼저 데이터프레임을 Spark SQL에 등록해야 합니다.
# Register df with Spark SQL
trios_with_referenceRegion.createOrReplaceTempView(“trios”)
이제 데이터프레임이 등록되었으므로 해당 데이터프레임에 대해 SQL 쿼리를 실행할 수 있습니다. 첫번째 쿼리의 경우 대체(ALT) 대립 유전자가 하나 이상 있는 샘플 NA 19685에 속하는 변형의 이름을 선택합니다.
# filter by alleles. This is a list of variant names that have an alternate allele for the child
alternate_variant_sites = spark.sql(“SELECT variant.names[0] AS snp FROM trios \
WHERE array_contains(alleles, ‘ALT’) AND sampleId == ‘NA19685′”)
collected_sites = map(lambda x: x.snp, alternate_variant_sites.collect())
다음 쿼리를 위해 부모가 둘 다 참고 대립 유전자를 가지고 있는 사이트를 필터링합니다. 그런 다음 이전에 아이로부터 생성된 세트를 기준으로 이러한 변형을 필터링합니다.
# get parent records and filter by only REF locations for variant names that were found in the child with an ALT
filtered1 = spark.sql(“SELECT * FROM trios WHERE sampleId == ‘NA19661’ or sampleId == ‘NA19660’ \
AND !array_contains(alleles, ‘ALT’)”)
filtered2 = filtered1.filter(filtered1[“variant.names”][0].isin(collected_sites))
snp_counts = filtered2.groupBy(“variant.names”).count().collect()
# collect snp names as a list
snp_names = map(lambda x: x.names, snp_counts)
denovo_snps = [item for sublist in snp_names for item in sublist]
denovo_snps[:10]
흥미로운 변형들을 발견했기 때문에, 메모리로부터 유전자형들을 유지시키지 않아도 됩니다.
trio_df.unpersist()
선형 데이터 작업
많은 잠재적인 새로운 변형 사이트들을 발견했습니다. 그런 다음 이러한 사이트 중 일부를 시각적으로 확인하여 원시 선형이 이러한 새로운 히트(hits)와 일치하는지 확인할 수 있습니다.
먼저, NA 19685 트리오에 대한 선형 데이터에 로드합니다.
# load in NA19685 exome from s3a
childReadsPath = ‘s3a://1000genomes/phase1/data/NA19685/exome_alignment/NA19685.mapped.illumina.mosaik.MXL.exome.20110411.bam’
parent1ReadsPath = ‘s3a://1000genomes/phase1/data/NA19685/exome_alignment/NA19660.mapped.illumina.mosaik.MXL.exome.20110411.bam’
parent2ReadsPath = ‘s3a://1000genomes/phase1/data/NA19685/exome_alignment/NA19661.mapped.illumina.mosaik.MXL.exome.20110411.bam’
childReads = ac.loadAlignments(childReadsPath, stringency=LENIENT)
parent1Reads = ac.loadAlignments(parent1ReadsPath, stringency=LENIENT)
parent2Reads = ac.loadAlignments(parent2ReadsPath, stringency=LENIENT)
이 예에서는 S3:// 스타일 URL 대신 S3a://를 사용합니다. 그 이유는 ADAM 형식이 Java NIO를 사용하여 BAM 파일에 액세스하기 때문입니다. 이를 위해 Hadoop 분배 파일 시스템에 대한 JSR 203 구현을 사용하여 이러한 파일에 액세스 합니다. 이것은 그 자체로 S3a:// 프로토콜을 필요로 합니다. 그 구현을 이 GitHub 저장소에서 볼 수 있습니다.
이제 트리오에서 각 3개의 개별 데이터 선형 데이터가 제공됩니다. 그러나 데이터는 아직 메모리에 로드되지 않았습니다. 데이터에 대한 빠른 후속 액세스를 위해 이러한 데이터셋을 캐시 하려면 캐시() 기능을 실행합니다.
# cache child RDD and count records
# takes about 2 minutes, on 4 c3.4xlarge worker nodes
childReads.cache()
# Count reads in the child
childReads.toDF().count()
# Output should be 95634679
선형 데이터의 품질 관리
게놈 선형 데이터의 품질을 시각적으로 다시 정의하기 위한 일반적인 분석은 커버리지 분포를 확인하는 것입니다. 커버리지 분포를 통해 전체 샘플에 걸쳐 있는 읽기 커버리지를 파악할 수 있습니다.
그런 다음, 염색체 17의 아이 선형 데이터에 대한 샘플 커버리지 분포를 생성합니다.
# Calculate read coverage
# Takes 2-3 minutes
childCoverage = childReads.transform(lambda x: x.filter(x.contigName == “17”)).toCoverage()
childCoverage.cache()
childCoverage.toDF().count()
# Output should be 51252612
이제 커버리지 데이터가 계산되고 캐시되었으므로, 염색체 17의 커버리지 분포를 계산하고 커버리지 분포를 표시합니다.
# Calculate coverage distribution
# You can check the progress in the SparkUI by navigating to
# <PUBLIC_MASTER_DNS>:8088 and clicking on the currently running Spark application.
cd = CoverageDistribution(sc, childCoverage)
x = cd.plot(normalize=True, cumulative=False, xScaleLog=True, labels=”NA19685″)
보고 있는 데이터가 외부 데이터이므로 표준 데이터처럼 보입니다. 따라서, 여러분은 낮은 커버리지에서 많은 수의 광경을 볼 수 있고 100개 이상의 읽을 수 있는 적은 수의 게놈 위치를 볼 수 있습니다. 이제 커버리지를 완료했으므로, 이러한 데이터셋을 다음 분석을 위해 메모리 공간을 확보할 필요가 없습니다.
childCoverage.unpersist()
발단자에서 미센스 변형체의 사이트 보기
선형 데이터 및 필터링 변형을 확인한 후에는 발단자에 잠재적인 미센스 돌연변이를 가진 4개의 유전자가 있으며, 이 중에 YBX2, ZNF286B, KSR1, 및 GNA13이 포함됩니다. 아이와 부모의 원시 읽기 내용을 필터링하고 표시하여 이러한 사이트를 시각적으로 확인할 수 있습니다.
먼저, 아이가 읽는 것을 보세요. GNA13 변형(63052580-63052581)의 위치로 확대하면 다음과 같이 이질적인 T에서 A로 호출하는 걸 볼 수 있습니다.
# missense variant at GNA13: 63052580-63052581 (SNP rs201316886)
# Takes about 2 minutes to collect data from workers
contig = “17”
start = 63052180
end = 63052981
genomicRDD.ViewAlignments(childReads, contig, start, end)
실제로 이 위치에 변형이 있는 것처럼 보입니다. 대체 대립유전자 A가 포함된 이질적인 SNP가 있을 수 있습니다. 부모 데이터를 확인하여 이 변형이 부모에 나타나지 않는지 확인하십시오.
# view missense variant at GNA13: 63052580-63052581 in parent 1
contig = “17”
start = 63052180
end = 63052981
genomicRDD.ViewAlignments(parent1Reads, contig, start, end)
이것은 이 변형이 실제로 부모들이 아닌, 단지 발단자에만 존재하는 필터를 확인시켜줍니다.
요약
요약하자면, 이 글은 Amazon EMR에서 ADAM 및 Mango를 설정하고 실행하는 방법을 보여 주었습니다. Amazon S3의 공용 데이터셋의 1000 유전체 데이터셋을 탐색하기 위해 쌍방향적인 노트북 환경에 있는 툴을 사용하는 방법에 대해 보여드렸습니다. 1000 유전체 데이터 품질을 검사하고, 게놈에서 흥미로운 변형을 확인하며, 원시 데이터의 시각화를 통해 결과를 검증하기 위해 툴을 사용했습니다.
Mango에 대한 자세한 내용은 Mango 사용 설명서를 참조하십시오.
원문 URL: https://aws.amazon.com/ko/blogs/big-data/exploratory-data-analysis-of-genomic-datasets-using-adam-and-mango-with-apache-spark-on-amazon-emr/
** 메가존 TechBlog는 AWS BLOG 영문 게재글중에서 한국 사용자들에게 유용한 정보 및 콘텐츠를 우선적으로 번역하여 내부 엔지니어 검수를 받아서, 정기적으로 게재하고 있습니다. 추가로 번역및 게재를 희망하는 글에 대해서 관리자에게 메일 또는 SNS페이지에 댓글을 남겨주시면, 우선적으로 번역해서 전달해드리도록 하겠습니다.