# 07. 한 샘플에서 아웃라이어 유전자 찾기

> 환자 한 명의 발현량 벡터를 GTEx·PCAWG 레퍼런스와 비교해 비정상적으로 튀는 유전자를 찾는 원리. 백분위·z-score·fold change·max 비교와 기능 분석의 해석 한계까지 다룬다.

RNA-seq로 발현량 행렬을 얻은 뒤, 환자 한 명의 발현량 벡터 하나만 있을 때, 공개 레퍼런스를 이용해 **비정상적으로 튀는 유전자를 어떻게 고르는지**에 집중합니다.

> **이 노트의 목표**: N=1 분석에서 비교 기준이 왜 필요한지, 백분위·z-score·fold change·max 비교가 각각 무엇을 말하는지, 그리고 "발현이 튄다"를 어디까지 해석할 수 있는지 이해한다.

운영 중인 서버 한 대의 지표만 보고 이상 여부를 판단할 수 없는 상황과 같습니다. CPU 80%가 높은지 알려면 같은 역할의 서버가 보통 어느 범위에서 움직이는지 알아야 합니다. 여기서는 서버 지표 벡터 대신 **샘플의 유전자 발현 벡터**, 과거 관측 분포 대신 **정상조직·암 레퍼런스 코호트**를 사용합니다.

```python
for gene in common_genes:
    baseline = reference_samples[gene]
    score[gene] = compare(patient_sample[gene], baseline)

candidates = rank(score)
```

이 의사코드의 출력은 원인 판정이 아니라 추가로 조사할 후보 순위입니다.

:::tip[실제 파일과 명령을 확인하고 싶다면]
STAR·RSEM·featureCounts 실행법과 로그·결과 파일 점검은 [RNA-seq 처음부터 끝까지 실습](/practice/rna-seq-from-fastq/)에서 다룹니다.
:::

## 1. 한 샘플로 무엇을 알 수 있나: 샘플이 딱 하나일 때

지금까지 만든 것은 **환자 한 명의 종양 조직 RNA-seq 프로파일 하나** 입니다. 유전자 수만 개(주석 기준 약 6만 개) 각각의 발현량(TPM 또는 raw count)이 담긴 **벡터 하나**.

DataFrame으로 보면 행은 유전자, 열은 샘플, 값은 발현량입니다. 환자 한 명은 열 하나이고, 같은 유전자의 레퍼런스 값들은 가로로 긴 분포를 만듭니다.

| gene | patient | normal_001 | normal_002 | cancer_001 | ... |
| --- | ---: | ---: | ---: | ---: | --- |
| `MDM2` | 1235.72 | 4.10 | 5.32 | 12.80 | ... |
| `ACTB` | 850.31 | 790.20 | 910.44 | 822.17 | ... |

> 그래서… **이 중 뭐가 중요한가?** 발현량 6만 줄을 눈으로 다 볼 수는 없습니다. 필요한 질문은 하나로 좁혀집니다:
>
> **"이 종양에서, 정상 조직이나 다른 암과 비교했을 때, 비정상적으로 발현이 튀는 유전자는 무엇인가?"**

이 질문에 답하는 작업이 **아웃라이어 유전자(outlier gene) 탐색**입니다.

![아웃라이어 개념: PCAWG 전체 암 분포에서 MDM2가 index tumor에서 크게 튀고, 대조로 ACTB는 분포 한가운데에 위치](/images/lessons/rna-seq-analysis/outlier-concept-mdm2.png)

*아이디어. 가로축은 한 유전자의 발현량(TPM), 세로축은 그 발현량을 가진 **샘플 수**. 대부분의 종양(All tumors)은 `MDM2`가 낮게 몰려 있는데, **이 환자(Index tumor)** 는 저 멀리 오른쪽: 육종(Sarcoma) 샘플들과 비교해도 튄다. 오른쪽 작은 그림의 `ACTB`(housekeeping 유전자)는 대조로, 누구에게나 높아 분포 한가운데 위치한다: 이런 건 아웃라이어가 아니다.*

### 튄다고 곧 원인은 아니다: 해석의 함정

발현량이 높게 튀는 걸 발견해도, **그 자체로 "이 유전자가 암을 일으켰다"고 말할 수는 없습니다.** 높은 발현은 다음 중 무엇이든 될 수 있습니다:

| 성격 | 의미 |
| --- | --- |
| **원인(cause)** | 종양 발생·진행에 **직접 기여** |
| **결과(consequence)** | 종양 관련 [pathway](/lessons/pathway/)가 켜지면서 **부차적으로** 함께 오른 발현 |
| **상태(state)** | 종양 미세환경 변화(저산소·영양 부족·전이 능력 획득 등)에 **적응**하며 바뀐 발현 |

> 지금 데이터에서 확실히 아는 건 **"발현이 튄다"** 뿐입니다. 인과를 밝히려면 아래가 **함께** 필요합니다:
>
> **WGS 변이 분석 + pathway 분석 + 문헌조사**: 예컨대 *과발현이 이 유전자 자체의 변이 때문인지*, 아니면 *변이는 없는데 다른 조절 유전자의 영향으로 튀는 것인지*를 가르고, 후자라면 그 조절 고리를 끊는 접근을 생각해 볼 수 있습니다. 여기에 **다른 오믹스**(예: 전사인자·효소가 붙을 수 있는 열린 자리를 보는 ATAC-seq)를 얹어 층층이 확인합니다.

그럼에도 아웃라이어 탐색이 유용한 이유는, **정상·타 환자 대비 변화가 분명한 유전자**를 추려 **표적 치료 후보·후보 약물 탐색**의 출발점으로 삼을 수 있기 때문입니다: Sid 케이스가 [MDM2 과발현](/start/case-sid/)을 표적으로 잡은 것처럼.

## 2. 비교 기준: 레퍼런스 데이터와 검증

한 샘플을 "비교"하려면 기준이 될 **공개 레퍼런스**가 필요합니다. Sid가 쓴 것과 같은, **바로 받아 쓸 수 있는 TPM 수준**의 두 데이터를 썼습니다:

| 레퍼런스 | 규모 | 성격 |
| --- | --- | --- |
| **GTEx** | 74,628 genes · 19,616 샘플 (31 조직) | **정상 조직** 발현 |
| **PCAWG** (EBI 경유) | 56,717 genes · 1,350 샘플 (27 암종) | **다른 암 환자** 발현 |

세 데이터셋(내 샘플 + GTEx + PCAWG)에 **모두 존재하는 공통 유전자 45,698개**만 추려서 비교했습니다.

이 교집합은 테이블 조인의 키를 맞추는 단계입니다. 데이터셋마다 유전자 ID 버전과 포함 범위가 다르므로 이름만 보고 열을 이어 붙이면 다른 유전자를 비교하거나 결측치를 대량으로 만들 수 있습니다. 유전자 ID 체계, 주석 버전, 발현량 단위와 전처리를 먼저 맞춰야 합니다.

> 🧩 그림의 암 레퍼런스는 **EBI(유럽 생물정보학 연구소)** 의 Expression Atlas에서 받은 PCAWG 데이터이며, 그림·표에는 **`ebi`** 라는 라벨로 표기했습니다. [Sid 사례](/start/case-sid/)에 나온 PCAWG와 같은 계열입니다.

### 검증 ①: Boston vs Tempus 상관

본론에 들어가기 전에, **만든 count가 믿을 만한지** 먼저 봅니다. 같은 검체를 두 회사에 맡긴 데이터이므로, 둘의 발현량은 **비슷한 분포**여야 합니다.

![Boston Gene(y) vs Tempus(x) 종양 TPM 산점도: 대각선을 따라 강한 상관, MDM2·PANX3는 오른쪽 위, FAP는 중간](/images/lessons/rna-seq-analysis/boston-tempus-validation.png)

*같은 환자·같은 시점(T0)을 두 회사에서 시퀀싱한 결과의 log2(TPM+1) 산점도. 완전한 직선은 아니지만 **대각선을 따라 강한 상관**을 보여, 발현량 정량이 재현성 있게 만들어졌음을 확인. 상관과 로그 변환의 의미는 [분포와 요약 통계](/reference/distributions-and-summary/)에서 설명한다.*

### 검증 ②: 이미 알려진 유전자가 실제로 튀는가

Sid 케이스에서 **질병 관련 유전자로 명시**됐던 것들이, 내 데이터에서도 같은 패턴으로 튀는지 봅니다.

![FAP 유전자의 TPM 분포: GTEx(정상, 파랑)와 EBI(암, 초록) 분포 위에 종양 값(주황 점선 Tumor=52.8)이 오른쪽에 위치](/images/lessons/rna-seq-analysis/fap-tpm-distribution.png)

*한 플롯이 **유전자 하나**에 대응한다. 가로축은 log2(TPM+1), 세로축은 **샘플 수**. 파란색은 GTEx(정상), 초록색은 EBI(암) 분포이고, **주황 점선이 내 종양 값**이다. `FAP`는 원래 싱글셀 데이터에서 찾아낸 표적이라 total RNA에서는 신호가 아주 도드라지진 않지만(그래도 과발현 쪽), 주요 표적으로 거론된 다른 유전자들은 분포를 한참 벗어나 찍힌다.*

> 읽는 법: 대부분의 샘플은 히스토그램 **중심** 근처에 발현량이 몰려 있고, **내 샘플(주황 선)이 그 중심에서 얼마나 오른쪽으로 벗어났는지**가 곧 "튀는 정도"입니다. 이건 Sid가 [공개 뷰어](https://osteosarc.com/rnaseq/)에서 보여준 그림과 같은 형식입니다.

## 3. 아웃라이어를 고르는 지표

유전자 4.5만 개를 눈으로 다 볼 수는 없으니, **수치로 아웃라이어를 추리는 지표**가 필요합니다. 종양 TPM을 레퍼런스 분포와 견주는 네 가지를 씁니다.

각 지표는 서로 다른 정렬 키입니다. 백분위는 분포 안의 위치, z-score는 평균으로부터의 표준화된 거리, fold change는 기준값과의 배수, max comparison은 관측 범위를 넘어섰는지를 나타냅니다. 한 열만 정답으로 두기보다 여러 열을 함께 보고 후보를 좁힙니다.

| 지표 | 무엇을 보나 | 기준 예시 |
| --- | --- | --- |
| **백분위 순위(percentile)** | 내 TPM이 레퍼런스 분포에서 **상위 몇 %** 인가 | max · 상위 1% · 5% 등 |
| **Z-score** | 레퍼런스 **평균에서 표준편차 몇 배** 떨어졌나 | **∣Z∣ ≥ 2** |
| **Fold change** | 레퍼런스 **중앙값(median)의 몇 배** 인가 | **log2FC ≥ 2** → 4배 이상 |
| **Max comparison** | 정상에서 관찰된 **최댓값보다도** 종양이 높은가 | True / False (직접 비교) |

Z-score는 레퍼런스의 평균에서 표준편차 몇 칸 떨어졌는지 나타냅니다. Fold change는 기준보다 몇 배 큰지 나타내며, log2FC 2는 4배, log2FC 9는 512배입니다. 백분위·평균·중앙값·표준편차·로그 변환의 계산과 한계는 [분포와 요약 통계](/reference/distributions-and-summary/)를 참고하세요.

Max comparison은 논문에 널리 통용되는 표준 지표라기보다 Sid 사례를 재현하기 위한 직접 비교입니다. “정상에서 관찰된 최댓값까지 넘었는가”를 True 또는 False로 기록합니다.

### 실제 계산: 알려진 표적 유전자들

암(EBI/PCAWG) 레퍼런스로 계산한 결과의 일부입니다. `ebi_median`·`ebi_max`는 암 환자 약 1,350명의 그 유전자 발현량의 **중앙값·최댓값**입니다.

| gene | tumor TPM | ebi median | ebi max | z-score | log2FC | ≈ 배수 |
| --- | --- | --- | --- | --- | --- | --- |
| **PANX3** | 619.43 | 0.138 | | **19.84** | 9.14 | **≈ 560×** |
| **MDM2** | 1235.72 | 4.954 | 2731 | 4.47 | 5.32 | ≈ 40× |
| **FAP** | 52.76 | 2.0 | 301 | 1.97 | 3.75 | ≈ 13× |

- `PANX3`를 보면: 암 환자 중앙값이 0.138인데 내 종양은 **619.43**. **백분위로 최댓값마저 넘고**, z-score는 **19.84**(평균에서 표준편차 20배 거리!), 중앙값 대비 **약 560배**(log2FC 9.14) 발현입니다. 교과서적인 아웃라이어죠.
- `MDM2`는 Sid가 실제 표적으로 삼은 유전자로, 종양 TPM 1235.72: 공개 뷰어에서는 1500 언저리로 표시되는 값입니다.
- `FAP`는 신호가 상대적으로 약합니다(z ≈ 2). total RNA에서는 도드라지지 않는, 싱글셀에서 잡힌 표적이라는 점과 일치합니다.

### 1차 필터: 정상/암 최댓값 대비 2배 선

여러 지표를 **조합**해 후보를 좁힙니다. 예컨대 종양 TPM을 **EBI 암 max TPM과 직접 산점도**로 놓고, 대각선(y=x)과 그 **2배 선**(log2 기준 +1)을 그으면:

![Boston 종양 TPM(y) vs EBI 암 max TPM(x) 산점도: y=x 대각선과 2배 선 위로 PANX3·MDM2가 크게 벗어나 있음](/images/lessons/rna-seq-analysis/max-comparison-scatter.png)

*점 하나가 유전자 하나(27개 암종·1,350 샘플의 max 기준). 대각선(y=x)은 "암 최댓값과 같음", 위쪽 점선은 그 **2배 선**. `PANX3`·`MDM2`처럼 **2배 선을 넘어 위로 벗어난 유전자**가 1차 아웃라이어 후보다. max 기준이라 종양에서만 튀는 완벽한 그림은 아니고(샘플 간 변이가 있어), 그래서 **여러 지표를 함께** 본다.*

> 이렇게 **레퍼런스 max·중앙값 대비 2배 이상 벗어난 유전자**를 1차 아웃라이어로 추려 냅니다. 한 지표로 순위를 매길 수도, 여러 지표를 조합할 수도 있습니다.

## 정리

한 샘플의 값은 혼자서는 높고 낮음을 말할 수 없습니다. 같은 유전자의 정상조직·다른 암 분포 위에 놓고 거리를 측정해야 합니다.

> **발현량 행렬 한 열 → 공통 유전자 선택 → 레퍼런스 분포와 비교 → 백분위·z-score·fold change로 후보 정렬**

아웃라이어는 추가로 조사할 후보이지 암의 원인이나 치료 표적이라는 결론이 아닙니다. 다음 [유전자 목록의 공통 기능 찾기](/lessons/gene-list-enrichment/)에서는 후보 유전자 목록에 어떤 기능이 몰려 있는지 해석합니다.

---

### 출처

- 정렬 도구: [STAR](https://github.com/alexdobin/STAR) · [Samtools](https://github.com/samtools/samtools)
- 정량: [RSEM](https://github.com/deweylab/RSEM) · [featureCounts (Subread)](https://subread.sourceforge.net/)
- 정상 조직 발현 레퍼런스: [GTEx Portal](https://gtexportal.org/home/downloads/adult-gtex/bulk_tissue_expression)
- 암 발현 레퍼런스: [PCAWG (*Nature* 2020)](https://www.nature.com/articles/s41586-020-1969-6) · [EBI Expression Atlas](https://www.ebi.ac.uk/gxa/home)
- 기능 강화 분석: [g:Profiler](https://biit.cs.ut.ee/gprofiler/gost) · [NetworkAnalyst](https://www.networkanalyst.ca/)
- Sid 공개 발현 뷰어: [osteosarc.com/rnaseq](https://osteosarc.com/rnaseq/)