LOCRC 실습 2: GO와 GSEA로 기능·경로 찾기
1. 무엇을 실습하나
섹션 제목: “1. 무엇을 실습하나”이 실습에서 답할 질문은 하나입니다.
수천 개의 차등 발현 결과에서 반복되는 생물학적 기능과 pathway를 어떻게 찾을까?
LOCRC 실습 1에서는 GSE251845의 Tumor와 인접 Normal 조직을 비교했습니다. 22명의 환자에게서 얻은 두 조직을 paired DESeq2로 분석한 결과, 저발현 필터를 통과한 유전자는 24,116개였습니다. 그중 padj < 0.05와 |log2FoldChange| >= 1을 함께 만족한 후보는 8,256개였습니다.
유전자 8,256개의 이름을 차례로 읽어서는 결과의 공통 구조를 파악하기 어렵습니다. 이 페이지에서는 유전자를 하나씩 보는 대신, 이미 알려진 기능별 유전자 집합과 대조합니다. 선택된 DEG 목록을 사용하는 over-representation analysis(ORA)와 전체 유전자 순위를 사용하는 gene set enrichment analysis(GSEA)를 차례로 실행합니다.
앞 실습에서 이어받는 입력
섹션 제목: “앞 실습에서 이어받는 입력”같은 colorectal-deg 프로젝트에서 계속 진행합니다. 다음 파일이 있어야 합니다.
| 입력 | 내용 | 이 페이지에서의 쓰임 |
|---|---|---|
outputs/GSE251845_tumor_vs_normal_PyDESeq2.csv | 24,116개 유전자의 log2FoldChange, padj 등 | DEG 선택, ORA background, GSEA 순위 |
outputs/ensembl-to-symbol.csv | Ensembl gene ID와 human gene symbol의 대응 | 분석 도구와 결과 그림의 ID 통일 |
analysis.py | 데이터 검증부터 DEG 후보 확인까지의 코드 | GO·GSEA 코드를 이어서 추가 |
파일이 없다면 LOCRC 실습 1의 프롬프트를 먼저 완료합니다. 이 페이지는 원본 count에서 DESeq2를 다시 설명하지 않고, 유전자별 결과표를 기능별 결과표로 바꾸는 과정에 집중합니다.
같은 결과표에서 두 가지 입력을 만든다
섹션 제목: “같은 결과표에서 두 가지 입력을 만든다”ORA와 GSEA는 같은 DESeq2 결과에서 출발하지만 서로 다른 모양의 입력을 받습니다.
| 방법 | 입력으로 만드는 것 | 한 줄로 묻는 질문 |
|---|---|---|
| ORA | 기준을 통과한 유전자 목록과 background | DEG 기준을 통과한 Tumor-up 목록에 세포 주기 관련 유전자가 예상보다 많이 들어 있는가? |
| GSEA | 모든 유전자와 순위 점수 | 예를 들어 세포 주기 관련 유전자들이 순위 위쪽이나 아래쪽에 함께 몰리는가? |
이 실습의 ORA 입력은 Tumor-up 3,631개와 Tumor-down 4,625개입니다. GSEA 입력은 경계를 통과했는지와 관계없이 순위를 만들 수 있는 전체 유전자입니다. 두 방법의 입력 차이를 이해해야 서로 다른 결과가 나온 이유를 설명할 수 있습니다.
경계 하나가 ORA와 GSEA의 입력을 나누는 장면
섹션 제목: “경계 하나가 ORA와 GSEA의 입력을 나누는 장면”아래는 방법의 차이를 보기 위한 가상 결과입니다.
| gene | log2FoldChange | padj | ORA의 Tumor-up 목록 | GSEA의 위치 |
|---|---|---|---|---|
GENE_A | 2.4 | 0.001 | 포함 | 위쪽 |
GENE_B | 1.1 | 0.020 | 포함 | 위쪽 |
GENE_C | 0.9 | 0.030 | 제외 | 위쪽에 가까움 |
GENE_D | 0.4 | 0.300 | 제외 | 중앙에 가까움 |
GENE_E | -1.3 | 0.010 | Tumor-down에 포함 | 아래쪽 |
GENE_B와 GENE_C의 변화량 차이는 작지만, |log2FoldChange| >= 1이라는 경계 때문에 ORA에서는 서로 다른 목록에 놓입니다. GSEA는 둘 다 순위 위쪽에 남깁니다. ORA는 기준을 통과한 DEG 목록을 설명하기 쉽고, GSEA는 경계 바로 밖에 있는 유전자까지 같은 방향의 신호에 포함할 수 있습니다.
ORA의 질문을 풀어서 보면
섹션 제목: “ORA의 질문을 풀어서 보면”ORA는 다음을 묻습니다.
예를 들어, 내가 정한 DEG 기준을 통과한 Tumor-up 유전자 목록에 세포 주기 관련 유전자가 무작위 선택에서 기대되는 수보다 많이 들어 있는가?
이 질문을 하는 이유는 유전자 하나의 결과보다 여러 유전자가 공유하는 기능이 반복되는지를 보기 위해서입니다.
100개의 공이 든 상자를 생각해 봅시다. 공 하나는 유전자 하나이고, 그중 10개에는 세포 주기라는 태그가 붙어 있습니다. 여기서 Tumor-up 후보 20개를 골랐습니다.
- 무작위로 20개를 골랐다면 세포 주기 공은 평균적으로 약 2개가 들어갑니다.
- 실제 후보 20개에 세포 주기 공이 8개 들어 있다면, 이 태그가 후보 목록에 유난히 많이 모인 셈입니다.
- 실제로 1개나 2개만 들어 있다면, 세포 주기 태그가 특별히 많이 모였다고 말하기 어렵습니다.
ORA는 이 차이가 우연히도 생길 만한 크기인지 통계적으로 계산합니다. 질문에 대한 답에 따라 먼저 내릴 수 있는 결론은 다음처럼 달라집니다.
| ORA 결과 | 먼저 말할 수 있는 결론 | 아직 말할 수 없는 것 |
|---|---|---|
| Tumor-up에서 세포 주기 term이 유의함 | Tumor에서 크게 증가한 후보 목록에 세포 주기 관련 유전자가 예상보다 많이 포함됨 | 세포 주기가 실제로 몇 배 빨라졌음 |
| Tumor-down에서 면역 term이 유의함 | Tumor에서 감소한 후보 목록에 면역 관련 유전자가 예상보다 많이 포함됨 | 암세포 내부의 면역 기능이 직접 억제됐음 |
| 해당 term이 유의하지 않음 | 현재 후보 기준으로 고른 목록에서는 그 기능이 과대표현됐다는 근거가 부족함 | 그 기능이 생물학적으로 전혀 관련 없음 |
마지막 행이 중요합니다. ORA에서 유의하지 않은 이유는 정말 공통 기능이 없어서일 수도 있지만, 관련 유전자들이 후보 경계 바로 밖에 있어서 목록에서 빠졌기 때문일 수도 있습니다. 이때 전체 순위를 보는 GSEA가 다른 답을 줄 수 있습니다.
GSEA의 질문을 풀어서 보면
섹션 제목: “GSEA의 질문을 풀어서 보면”GSEA는 다음을 묻습니다.
DEG 기준을 통과한 유전자만 고르지 않고, 모든 유전자를 Tumor에서 가장 많이 증가한 유전자부터 가장 많이 감소한 유전자까지 순서대로 나열했을 때, 예를 들어 세포 주기 관련 유전자들이 순위의 앞쪽이나 뒤쪽에 함께 모여 있는가?
이 질문을 하는 이유는 변화가 아주 크지는 않아도 같은 방향으로 움직이는 유전자가 여러 개 있을 수 있기 때문입니다.
이번에는 100개 유전자를 Tumor에서 많이 증가한 순서부터 많이 감소한 순서까지 한 줄로 세웠다고 생각해 봅시다. 세포 주기 유전자 10개의 위치를 표시합니다.
- 10개 중 7개가 앞쪽 20칸에 모이면, 세포 주기 유전자들이 Tumor-up 방향으로 함께 움직이는 신호입니다.
- 10개 중 7개가 뒤쪽 20칸에 모이면, Normal에서 높은 방향으로 함께 움직이는 신호입니다.
- 10개가 순위 전체에 고르게 흩어지면, 한쪽 방향으로 모였다는 근거가 약합니다.
GSEA는 이렇게 한쪽으로 몰린 정도를 NES(Normalized Enrichment Score, 정규화된 enrichment 점수)로 요약합니다. NES는 분석 전체에 하나만 나오는 값이 아니라, 검사한 유전자 집합마다 하나씩 계산됩니다. 예를 들어 E2F targets, G2M checkpoint, Inflammatory response는 각각 별도의 NES를 가집니다. 유전자 집합마다 포함된 유전자 수가 달라 원래 점수인 ES를 그대로 비교하기 어렵기 때문에, 각 집합에서 무작위로 기대되는 점수 크기를 기준으로 보정한 값입니다. 이 페이지처럼 Tumor-up 유전자를 순위 앞쪽에 놓았다면, 양의 NES는 Tumor-up 쪽에 몰렸다는 뜻이고 음의 NES는 Tumor-down 쪽에 몰렸다는 뜻입니다. 절댓값이 클수록 한쪽으로 몰리는 경향이 강하지만, 통계적으로 유의한지는 NES만 보지 않고 FDR q-value도 함께 확인합니다.
따라서 곡선의 피크가 곧 NES인 것은 아닙니다. 누적 enrichment score 곡선이 0에서 가장 멀어진 최고점 또는 최저점이 ES이고, 이 ES를 무작위 결과와 비교해 정규화한 값이 NES입니다.
이 실습은 큰 양의 log2FoldChange부터 큰 음의 값까지 내림차순으로 정렬합니다. 따라서 GSEA의 답은 다음처럼 읽습니다.
| GSEA 결과 | 먼저 말할 수 있는 결론 | 아직 말할 수 없는 것 |
|---|---|---|
| 유의한 양의 NES | gene set 구성원들이 Tumor-up 쪽에 집단적으로 몰림 | 해당 pathway가 실제로 활성화됐음 |
| 유의한 음의 NES | gene set 구성원들이 Normal-up, 즉 Tumor-down 쪽에 집단적으로 몰림 | 해당 pathway가 암세포 안에서 직접 억제됐음 |
| NES가 0에 가깝거나 유의하지 않음 | 현재 순위에서는 구성원들이 한쪽에 모였다는 근거가 부족함 | 구성원 중 중요한 개별 유전자가 하나도 없음 |
GSEA의 양수와 음수는 고정된 생물학적 뜻이 아닙니다. 어떤 점수로 어느 방향부터 정렬했는지가 부호의 뜻을 결정합니다. 이 페이지에서는 Tumor-up부터 정렬했기 때문에 양의 NES가 Tumor-up을 가리킵니다.
두 질문의 답을 함께 보면
섹션 제목: “두 질문의 답을 함께 보면”ORA와 GSEA는 같은 답을 반복하는 검사가 아닙니다. 하나는 선택된 후보 목록의 비율, 다른 하나는 전체 순위에서의 위치를 봅니다.
| ORA | GSEA | 먼저 확인할 해석 |
|---|---|---|
| 유의함 | 같은 방향으로 유의함 | 큰 변화를 통과한 후보와 전체 순위가 같은 기능 신호를 함께 지지함 |
| 유의하지 않음 | 유의함 | 큰 변화의 후보는 적지만, 작은 변화의 여러 유전자가 같은 방향에 넓게 모였을 수 있음 |
| 유의함 | 유의하지 않음 | 경계를 통과한 일부 유전자의 겹침은 크지만 gene set 전체가 한 방향으로 움직이지 않을 수 있음 |
| 유의하지 않음 | 유의하지 않음 | 현재 후보 기준과 순위에서는 그 gene set의 집단적 변화 근거가 부족함 |
예를 들어 세포 주기가 Tumor-up ORA와 양의 NES GSEA에서 모두 나타나면, “크게 증가한 후보에도 세포 주기 유전자가 많고, 경계 밖 유전자까지 포함한 전체 순위에서도 Tumor-up 쪽에 몰린다”고 정리할 수 있습니다. 이것이 두 분석이 같은 방향을 가리킨다는 뜻입니다. 그래도 RNA 발현 결과만으로 단백질 활성이나 실제 세포 분열 속도를 확정하지는 않습니다.
2. 유전자를 기능별 집합으로 바꾸어 보기
섹션 제목: “2. 유전자를 기능별 집합으로 바꾸어 보기”DESeq2 결과표의 한 행은 유전자 하나를 설명합니다. 설명을 위해 CDK1, CDC20, MCM2가 모두 Tumor에서 증가한 세 행이 있다고 가정해 봅시다. 이 결과만 보고는 세 유전자가 어떤 일을 함께 하는지 바로 알기 어렵습니다. 그래서 각 유전자 ID를 이미 만들어진 기능 목록과 대조합니다.
| 설명용 유전자별 결과 | 기능 자료에서 찾은 목록의 예 |
|---|---|
CDK1 증가 | GO의 cell cycle, Reactome의 Cell Cycle, Hallmark의 E2F_TARGETS |
CDC20 증가 | GO의 cell cycle, Reactome의 Cell Cycle, Hallmark의 E2F_TARGETS |
MCM2 증가 | DNA 복제 관련 GO term과 Reactome pathway, Hallmark의 E2F_TARGETS |
이 과정은 유전자 자체를 다른 것으로 바꾸는 것이 아닙니다. 유전자별 결과표를 기능 자료의 소속표와 대조한 뒤, 같은 항목에 연결된 유전자들을 다시 모으는 과정입니다. 한 유전자는 여러 기능에 참여할 수 있으므로 여러 목록에 동시에 들어갈 수 있습니다.
Gene Ontology: 유전자에 붙은 기능 태그를 모은다
섹션 제목: “Gene Ontology: 유전자에 붙은 기능 태그를 모은다”Gene Ontology의 GO:BP는 유전자 산물이 어떤 생물학적 과정에 참여하는지 표준 용어로 기록합니다. 예를 들어 GO:0007049 cell cycle은 세포가 유전물질을 복제하고 나누는 전체 과정에 붙은 이름입니다.
CDK1, CCNB1, CDC20 같은 유전자가 이 term과 연결돼 있다고 해봅시다. Tumor-up 목록에 이 유전자들이 많이 들어 있으면 ORA는 “Tumor-up 유전자에 cell cycle 태그가 예상보다 자주 나타나는가?”를 검사합니다.
GO는 같은 일에 참여한다는 분류에 강합니다. 다만 CDK1 다음에 CDC20이 작동한다는 식의 반응 순서를 보여 주는 지도는 아닙니다. 또한 cell cycle 아래에 mitotic cell cycle, chromosome segregation 같은 더 구체적인 term이 연결되므로, 비슷한 결과가 여러 줄 함께 나올 수 있습니다.
Reactome: 연결된 반응을 하나의 pathway로 묶는다
섹션 제목: “Reactome: 연결된 반응을 하나의 pathway로 묶는다”Reactome의 R-HSA-1640170 Cell Cycle은 세포 주기에 관여하는 분자와 반응을 연결한 pathway입니다. 이 지도에는 CCNB1과 CDK1이 복합체를 이루는 상태, CDC20이 관여하는 유사분열 단백질 분해 같은 구체적인 장면이 들어 있습니다.
Reactome을 ORA에 넣을 때는 이 복잡한 지도를 다음처럼 단순한 유전자 목록으로 펼쳐서 사용합니다.
R-HSA-1640170 Cell Cycle→ CDK1, CCNB1, CDC20, ...Tumor-up 목록과 이 목록이 많이 겹치면 Cell Cycle pathway가 유의하게 나올 수 있습니다. 이때 ORA가 사용한 정보는 어떤 유전자가 pathway에 속하는가입니다. 반응의 순서, 단백질의 변형 상태와 세포 내 위치까지 RNA-seq 결과가 확인했다는 뜻은 아닙니다.
MSigDB Hallmark: 함께 움직이는 대표 발현 신호를 모은다
섹션 제목: “MSigDB Hallmark: 함께 움직이는 대표 발현 신호를 모은다”MSigDB Hallmark는 여러 연구의 유전자 집합에서 반복되는 발현 신호를 대표 목록으로 정리한 컬렉션입니다. Reactome처럼 반응의 순서를 그린 지도가 아니라, 특정 생물학적 상태를 잘 나타내는 유전자 목록에 가깝습니다.
예를 들어 HALLMARK_E2F_TARGETS에는 E2F 전사인자의 표적으로 알려진 세포 주기 관련 유전자들이 들어 있습니다.
HALLMARK_E2F_TARGETS→ CDK1, CDC20, MCM2, PCNA, ...GSEA는 이 목록을 한 번에 하나씩 꺼내 전체 유전자 순위와 대조합니다. CDK1, CDC20, MCM2, PCNA 같은 구성원이 Tumor-up 쪽에 함께 몰리면 HALLMARK_E2F_TARGETS에 양의 NES가 계산될 수 있습니다. 다음에는 HALLMARK_G2M_CHECKPOINT를 꺼내 같은 계산을 처음부터 따로 수행합니다.
세 자료를 나란히 보면
섹션 제목: “세 자료를 나란히 보면”같은 세포 주기 신호를 다루더라도 세 자료가 붙이는 이름과 범위는 다릅니다.
| 자료 | 원래 묻는 관점 | 계산할 때 주로 사용하는 형태 | 이 예시에서 읽는 말 |
|---|---|---|---|
| GO:BP | 이 유전자는 어떤 과정에 참여하는가? | term마다 연결된 유전자 목록 | Tumor-up 목록에 cell cycle 관련 유전자가 많다 |
| Reactome | 어떤 분자 반응들이 하나의 흐름을 이루는가? | pathway에 포함된 유전자 목록 | Tumor-up 목록에 Cell Cycle pathway 유전자가 많다 |
| Hallmark | 어떤 대표 발현 신호의 유전자들이 함께 움직이는가? | 중복을 줄여 정리한 gene set | E2F_TARGETS 구성원이 Tumor-up 순위에 몰린다 |
따라서 세 자료에서 세포 주기 관련 결과가 모두 나와도 독립적인 발견 세 개로 바로 세지 않습니다. 서로 겹치는 유전자들이 같은 변화를 다른 지식 체계에서 반복해 보여 주는 것인지 실제 구성원 목록을 확인해야 합니다.
3. 선택한 DEG 목록으로 ORA 실행하기
섹션 제목: “3. 선택한 DEG 목록으로 ORA 실행하기”ORA는 선택한 목록과 하나의 기능 집합이 얼마나 많이 겹치는지 검사합니다. 한 term을 검사할 때 필요한 수는 네 가지입니다.
| 수 | 이 실습에서 가리키는 것 |
|---|---|
| background 크기 | DESeq2가 실제로 검정한 24,116개 유전자 |
| 후보 목록 크기 | Tumor-up 3,631개 또는 Tumor-down 4,625개 |
| term 크기 | background 중 특정 GO term이나 Reactome pathway에 연결된 유전자 수 |
| 교집합 크기 | 후보 목록과 term이 실제로 공유한 유전자 수 |
예를 들어 설명용 term이 background에서 200개 유전자와 연결돼 있다고 가정해 봅시다. Tumor-up 3,631개를 무작위로 골랐다면 단순 비율로 기대되는 교집합은 3,631 × 200 ÷ 24,116, 약 30개입니다. 실제 교집합이 80개라면 이 기능의 유전자가 후보 목록에 예상보다 많이 모였는지 통계적으로 검사할 이유가 생깁니다.
g:Profiler의 ORA는 이런 겹침을 누적 초기하분포 검정으로 계산하고, 수많은 term을 동시에 검사한 결과를 보정합니다. 이 실습에서 결과 열 이름은 p_value이지만 g:Profiler의 기본 g:SCS 다중검정 보정이 적용된 값입니다.
background는 가능한 후보의 범위다
섹션 제목: “background는 가능한 후보의 범위다”background는 단순히 알려진 모든 인간 유전자가 아닙니다. 이 실험에서 DEG가 될 기회가 있었던 유전자의 범위입니다. 저발현 필터에서 제외된 유전자는 후보 목록에 들어갈 수 없었으므로, 모든 인간 유전자를 background로 쓰면 후보가 될 수 없던 유전자까지 비교 대상에 섞입니다.
이 페이지에서는 DESeq2가 실제로 검정한 24,116개를 custom background로 사용합니다. g:Profiler에는 domain_scope="custom_annotated"와 background ID 목록을 함께 전달합니다. 이 설정은 custom background 안에서 해당 자료원에 annotation이 있는 유전자만 통계 범위로 사용합니다.
증가와 감소 목록을 나누는 이유
섹션 제목: “증가와 감소 목록을 나누는 이유”Tumor-up과 Tumor-down을 한 목록으로 합치면 “발현이 달라진 유전자에 면역 관련 기능이 많다”는 말은 할 수 있어도, 그 기능이 어느 방향의 유전자에서 나온 것인지 잃게 됩니다. 이 실습은 두 목록을 따로 분석해 다음 질문을 구분합니다.
- Tumor-up ORA: Tumor에서 증가한 유전자에 어떤 기능이 많이 포함됐는가?
- Tumor-down ORA: Tumor에서 감소한 유전자에 어떤 기능이 많이 포함됐는가?
1단계 프롬프트: g:Profiler GO·Reactome ORA
LOCRC 실습 1의 분석을 유지하고 Tumor-up과 Tumor-down DEG의 over-representation analysis를 추가해줘.
1. 현재 폴더가 colorectal-deg 프로젝트인지 확인하고, outputs/GSE251845_tumor_vs_normal_PyDESeq2.csv와 outputs/ensembl-to-symbol.csv가 있는지 검증해. 없으면 임의로 다시 만들지 말고 LOCRC 실습 1에서 빠진 파일을 보고하고 멈춰.
2. padj < 0.05이면서 abs(log2FoldChange) >= 1인 유전자만 사용해. 양수는 Tumor-up, 음수는 Tumor-down 목록으로 나눠.
3. 저발현 필터를 통과해 DESeq2에서 실제로 검정한 24,116개 유전자 전체를 ORA의 background로 사용해. g:Profiler 호출에 domain_scope="custom_annotated"와 background 목록을 명시하고 배경 유전자 수도 출력해.
4. gprofiler-official의 GProfiler를 사용하고 organism="hsapiens", sources=["GO:BP", "REAC"]로 각 목록을 따로 분석해.
5. 결과의 p_value에는 g:Profiler의 다중검정 보정이 적용됐다는 점을 코드 주석과 출력 설명에 명시해.
6. 각 결과에서 source, native, name, term_size, intersection_size, p_value, intersections 열을 outputs/ora-up.csv와 outputs/ora-down.csv에 저장해. intersections의 Ensembl ID는 그대로 보존하고, 앞 실습의 매핑을 이용한 intersection_symbols 열도 추가해.
7. 각 방향에서 p_value가 가장 작은 15개 term의 -log10(p_value)를 가로 막대그래프로 그려. GO:BP와 REAC의 색을 구분해 outputs/ora-up.png와 outputs/ora-down.png에 저장해.
8. 상위 term이 같은 유전자를 반복해서 공유하는지 intersections를 비교하고, 의미가 겹치는 term을 묶어 설명해.
9. uv run python analysis.py로 전체 스크립트를 실행하고 두 그래프를 직접 확인해. Tumor-up과 Tumor-down에서 각각 두드러지는 기능을 한국어로 설명한 뒤 멈춰.
유전자 목록의 기준을 바꾸면 ORA 결과도 바뀐다는 한계를 반드시 설명해.실제 실행에서는 Tumor-up 3,631개와 Tumor-down 4,625개를 따로 분석했습니다. 저발현 필터를 통과해 DESeq2가 검정한 24,116개 유전자를 background로 사용했으며, g:SCS 보정 뒤 유의한 term은 각각 197개와 395개였습니다. 결과표는 outputs/ora-up.csv와 outputs/ora-down.csv에 저장됐습니다.
현재 실험의 outputs/ora-up.png. Tumor에서 증가한 DEG에는 세포 주기, 염색체 분리와 DNA 복제 관련 유전자가 많이 포함됐습니다.
현재 실험의 outputs/ora-down.png. Tumor에서 감소한 DEG에는 다세포 생물 과정, 자극 반응, 세포 간 신호와 면역 관련 유전자가 많이 포함됐습니다.
막대가 길다는 것은 보정 p-value가 작다는 뜻이지, 발현량이 더 많이 변했다거나 과정이 더 강하게 작동했다는 뜻이 아닙니다. 변화량은 원래 DEG 표에서 확인하고, ORA 막대는 목록 안에서 기능 annotation이 얼마나 예상 밖으로 집중됐는지를 읽습니다.
비슷한 term 여러 개를 하나의 신호로 묶기
섹션 제목: “비슷한 term 여러 개를 하나의 신호로 묶기”Cell Cycle, Cell Cycle, Mitotic, Chromosome Segregation처럼 비슷한 이름이 함께 나타날 수 있습니다. GO는 상위·하위 term이 연결된 그래프이고 Reactome도 큰 pathway가 세부 pathway를 포함하므로, 이 term들은 같은 유전자 일부를 공유합니다.
결과표의 intersection_symbols를 비교해 같은 유전자가 막대 여러 개를 만들었는지 확인합니다. 겹침이 크다면 “세포 주기 계열의 신호”로 먼저 묶고, 그 안에서 현재 데이터가 뒷받침하는 가장 구체적인 term을 살펴봅니다. term 열다섯 개를 독립적인 발견 열다섯 개로 세지 않습니다.
4. 전체 유전자 순위로 GSEA 실행하기
섹션 제목: “4. 전체 유전자 순위로 GSEA 실행하기”ORA는 후보 기준이 필요합니다. |log2FoldChange| 경계를 1에서 0.8로 바꾸거나 padj 기준을 바꾸면 입력 목록도 달라지고 결과도 달라질 수 있습니다. GSEA는 이 문제를 다른 방식으로 봅니다. 유전자를 선택해 자르지 않고, 전체를 Tumor에서 높은 쪽부터 Normal에서 높은 쪽까지 정렬합니다.
이 실습에서는 log2FoldChange를 순위 점수로 사용합니다.
Tumor에서 높음 Normal에서 높음큰 양의 log2FC → 작은 양의 값 → 0 → 작은 음의 값 → 큰 음의 값한 gene set의 유전자가 왼쪽에 몰리면 Tumor-up 방향, 오른쪽에 몰리면 Normal-up 방향의 신호입니다. 순위 점수와 정렬 방향을 바꾸면 부호의 뜻도 바뀔 수 있으므로, 결과를 저장할 때 무엇으로 내림차순 정렬했는지 기록합니다.
running enrichment score는 순위표를 걷는 값이다
섹션 제목: “running enrichment score는 순위표를 걷는 값이다”GSEA는 순위표 맨 위에서 아래로 한 행씩 이동합니다.
- 현재 유전자가 검사할 gene set에 있으면 점수를 올립니다.
- gene set에 없으면 점수를 조금 내립니다.
- 순위의 양끝 중 한쪽에 구성원이 몰릴수록 0에서 멀리 벗어납니다.
- 0에서 가장 멀어진 지점이 enrichment score(ES)가 됩니다.
상위의 변화량이 큰 유전자는 더 큰 가중치를 받을 수 있습니다. 따라서 단순히 구성원 수만 세는 것이 아니라, gene set 구성원이 순위 어디에 있고 그 위치의 점수가 얼마나 큰지도 반영합니다.
ES, NES와 leading edge
섹션 제목: “ES, NES와 leading edge”| 결과 | 먼저 읽을 뜻 |
|---|---|
ES | 해당 gene set이 현재 순위의 한쪽에 몰린 정도 |
NES | gene set 크기 등의 영향을 보정해 비교하기 쉽게 만든 ES |
NOM p-val | permutation에서 현재 ES만큼 극단적인 값이 나온 비율 |
FDR q-val | 여러 gene set을 검사한 뒤의 false discovery rate 추정치 |
Lead_genes | ES가 peak에 도달하는 데 가장 크게 기여한 구성원 |
NES가 양수면 이 페이지의 내림차순 순위에서 Tumor-up 쪽에, 음수면 Normal-up 쪽에 gene set 구성원이 몰렸다는 뜻입니다. 이는 pathway가 켜지거나 꺼진 정도를 나타내는 단위가 아닙니다.
NOM p-val이 0으로 출력되더라도 확률이 정확히 0이라는 뜻은 아닙니다. 1,000번 permutation에서 더 극단적인 값을 한 번도 관측하지 못했다는 뜻이므로, 해상도는 permutation 횟수에 제한됩니다.
leading edge는 양의 ES에서는 peak까지 나타난 gene set 구성원, 음의 ES에서는 음의 peak를 만드는 반대쪽 구성원입니다. pathway 이름만 읽지 않고 실제 신호를 만든 유전자를 되짚을 때 사용합니다.
2단계 프롬프트: Hallmark GSEA와 running enrichment plot
이전 ORA 분석을 유지하고 전체 유전자 순위를 사용하는 Hallmark preranked GSEA를 추가해줘.
1. gseapy가 설치돼 있지 않으면 uv add gseapy를 실행하고, 설치된 버전을 출력해.
2. 앞 실습에서 저장한 outputs/ensembl-to-symbol.csv를 재사용해 Ensembl ID를 human gene symbol로 바꿔. 파일이 없을 때만 g:Profiler identifier conversion을 다시 실행해. 변환되지 않은 ID는 제외하고, 같은 symbol이 여러 번 나오면 abs(log2FoldChange)가 가장 큰 행 하나만 남겨.
3. 모든 유전자를 log2FoldChange 내림차순으로 정렬한 2열 rank table을 outputs/gsea-rank.tsv로 저장해. 첫 행과 마지막 행을 출력해 순위 방향을 검증해.
4. GSEApy prerank와 MSigDB_Hallmark_2020 gene set을 사용해. min_size=15, max_size=500, permutation_num=1000, threads=1, seed를 고정해 재현 가능하게 실행해.
5. 결과의 Term, ES, NES, NOM p-val, FDR q-val, Tag %, Gene %, Lead_genes를 outputs/gsea-hallmark.csv로 저장해.
6. FDR q-val이 작은 pathway를 골라 x축 NES의 dot plot을 outputs/gsea-hallmark-dotplot.png에 저장해. NES가 양수면 빨간색 계열, 음수면 파란색 계열로 표시하고 점 크기는 gene set size를 반영해.
7. E2F_TARGETS의 running enrichment score, gene hit 위치, ranked metric을 한 그림에 표시해 outputs/gsea-e2f-running-enrichment.png에 저장해.
8. 양의 NES 상위와 음의 NES 상위 pathway를 각각 출력하고, leading-edge gene도 함께 보여 줘.
9. uv run python analysis.py로 실행하고 두 그래프를 직접 확인해. ORA와 다른 점, NES 부호, E2F 곡선의 peak를 한국어로 설명한 뒤 멈춰.
온라인 gene set이나 identifier 변환에 접근하지 못하면 임의의 결과를 만들지 말고 어떤 단계에서 막혔는지 보고해.
GSE251845 Hallmark GSEA 결과표 예시.
GSE251845 Hallmark GSEA dot plot 예시.
오른쪽의 양의 NES에는 E2F targets, MYC targets, G2M checkpoint처럼 세포 주기와 증식에 관련된 gene set이 있습니다. 왼쪽의 음의 NES에는 fatty acid metabolism, interferon gamma response처럼 대사와 면역 반응에 관련된 gene set이 있습니다. 점의 색은 통계적 근거, 크기는 gene set의 크기를 나타내므로 x축 위치만 보지 않고 함께 읽습니다.
GSE251845 E2F_TARGETS running enrichment plot 예시.
가로축은 Tumor에서 높은 유전자부터 Normal에서 높은 유전자까지 정렬한 전체 목록입니다. 아래의 검은 세로선은 E2F_TARGETS 유전자가 순위에서 나타나는 위치입니다. 초반에 선이 몰리면서 초록색 running score가 빠르게 상승하므로, E2F 표적 유전자들이 Tumor-up 쪽에 집중됐다고 읽습니다.
peak 이전의 E2F 구성원은 양의 enrichment를 주로 만든 leading edge입니다. peak 이후에는 gene set에 속하지 않는 유전자를 지나면서 점수가 내려가고, 목록 끝에서는 다시 0에 가까워집니다.
5. ORA와 GSEA를 함께 읽기
섹션 제목: “5. ORA와 GSEA를 함께 읽기”두 방법은 우열 관계가 아니라 같은 결과표를 다른 질문으로 읽는 방법입니다.
| 확인할 점 | ORA | GSEA |
|---|---|---|
| 사용하는 유전자 | 기준을 통과한 후보 | 순위를 만들 수 있는 전체 유전자 |
| 변화 방향 | up과 down 목록을 따로 분석 | 순위 방향과 NES 부호로 표현 |
| 경계의 영향 | 큼 | 후보 경계는 사용하지 않음 |
| 대표적인 확인 대상 | intersections | Lead_genes와 hit 위치 |
| 답하는 질문 | DEG 기준을 통과한 목록에, 예를 들어 세포 주기 관련 유전자가 기대보다 많이 포함됐는가? | 예를 들어 세포 주기 관련 유전자들이 순위 한쪽에 몰리는가? |
현재 결과에서 Tumor-up ORA에는 세포 주기, 염색체 분리와 DNA 복제가 나타났고, GSEA의 양의 NES에도 E2F targets, MYC targets, G2M checkpoint가 나타났습니다. 서로 다른 입력 규칙에서 증식 관련 유전자 집합이 같은 Tumor-up 방향을 보였다는 점은 결과를 설명하는 중요한 공통 패턴입니다.
Tumor-down ORA와 음의 NES에는 면역 반응과 대사 관련 집합이 나타났습니다. 이 결과만으로 암세포 안에서 해당 pathway가 억제됐다고 결론 내릴 수는 없습니다. GSE251845는 여러 세포가 섞인 bulk 조직이므로, Tumor와 인접 Normal에 포함된 면역세포나 상피세포 비율의 차이도 같은 발현 패턴을 만들 수 있습니다.
ORA와 GSEA가 다르게 나오는 경우도 실패는 아닙니다. 먼저 다음을 확인합니다.
- ORA 후보 경계 바로 밖에 같은 방향의 유전자가 많이 있는가?
- 큰 변화의 소수 유전자가 아니라 작은 변화의 여러 유전자가 순위 한쪽에 모였는가?
- GO·Reactome·Hallmark가 같은 생물학적 개념을 서로 다른 구성원으로 정의했는가?
- ID 변환에서 제외되거나 하나의 symbol로 합쳐진 유전자가 많은가?
- 비슷한 term과 gene set이 같은 구성원을 반복해서 공유하는가?
결과를 재현하려면 입력 DEG 기준, background, 순위 점수와 방향, ID 변환 규칙, gene set 컬렉션 이름, 분석 날짜, g:Profiler와 GSEApy 버전, permutation 횟수와 seed를 함께 기록합니다.
6. 결과를 본 뒤 답해 볼 질문
섹션 제목: “6. 결과를 본 뒤 답해 볼 질문”- ORA background에 모든 인간 유전자 대신 실제 검정한 24,116개를 사용한 이유는 무엇인가?
- Tumor-up과 Tumor-down을 합치지 않고 따로 분석한 이유는 무엇인가?
- 비슷한 GO term 여러 개를 독립적인 발견으로 세면 안 되는 이유는 무엇인가?
- ORA에서 제외된 유전자도 GSEA 신호에 기여할 수 있는 이유는 무엇인가?
- 이 순위에서 양의 NES와 음의 NES는 각각 어느 조건을 가리키는가?
- running enrichment plot의 peak와 leading edge는 무엇을 보여 주는가?
- 면역 관련 gene set의 음의 NES를 암세포 안의 면역 억제로 바로 해석할 수 없는 이유는 무엇인가?
다음 조기·후기 발병 대장암 비교 실습에서는 같은 paired 분석을 두 연령군에서 따로 실행한 뒤, 유전자별 변화량을 직접 비교합니다.
공식 자료
섹션 제목: “공식 자료”- Gene Ontology 문서: GO term, 세 가지 분류와 관계 구조
- g:Profiler 문서: ORA, custom background와 다중검정 보정
- Reactome Pathway Analysis: pathway gene set과 over-representation analysis
- GSEA User Guide: ES, NES, FDR와 leading edge
- GSEApy 문서: Python
prerank실행 방법