Thuật toán PCA

PCA (Principal Component Analysis) hay Phân tích Thành phần Chính là một kỹ thuật giảm chiều dữ liệu, giúp tóm tắt và mô tả các đặc trưng chính của dữ liệu mà không mất nhiều thông tin quan trọng. Mục tiêu của PCA là tìm ra những hướng (gọi là thành phần chính) mà dữ liệu phân bố rộng nhất, tức là có phương sai lớn nhất.

Chuẩn hóa dữ liệu

Trừ đi giá trị trung bình và chia cho độ lệch chuẩn để khiến mọi chiều trong dữ liệu về cùng 1 thang đo * \(X centered=X−mean(X)\) * \(Xscaled = Xcentered/std(X)\)

Tìm giá trị trung tâm của dữ liệu

Gióng tất cả dữ liệu lên các trục và tính giá trị trung bình của từng trục.
Giao nhau của các đường giá trị trung bình trong không gian tại 1 điểm. Điểm đó là điểm chính giữa của dữ liệu

Di chuyển hệ trục

Di chuyển tập dữ liệu sao cho điểm chính giữa của dữ liệu trở thành gốc tọa độ

PCA1

Vẽ một đường thẳng bất kì đi qua gốc tọa độ mới
Xoay đường này sao cho khoảng cách của đường với các điểm dữ liệu gần nhất. Điều này đồng nghĩa với tổng khoảng cách từ gốc tọa độ đến hình chiếu của điểm dữ liệu lên trục là cao nhất (thông qua định lý Pythago)

PCA2

Tiếp tục áp dụng cách vẽ của PCA1, có được PCA2
Lấy 2 PCA này để làm trục tọa độ, chấp nhận mất các thông tin khác

Vẽ lại đồ thị

Chiếu các điểm dữ liệu lên 2 trục PCA1 và 2 Xoay hệ trục sao cho PCA1, tương ứng với trục x và PCA2 tương ứng với trục y.
Dựa vào các điểm chiếu trên 2 trục PCA để vẽ lại các điểm dữ liệu

Eigen value

Giá trị đo độ phân tán dữ liệu (dựa trên phương sai). Nếu giá trị này càng lớn thì PCA đó càng mang nhiều thông tin của dữ liệu
Khoảng cách từ gốc tọa độ đến hình chiếu của điểm lên PCA gọi là d
eigen value = tổng bình phương d/n-1 = variance
Như vậy eigen value = sd^2
Để tính khả năng đóng góp của PCR: \(sd^2/sum(sd)*100\)

Để hiểu rõ cách hoạt động của thuật toán, xem video sau Thuật toán PCA{https://www.youtube.com/watch?v=FgakZw6K1QQ&t=678s}

PCA in R

Tạo dữ liệu biểu hiện gen

# Tạo dữ liệu mẫu cho biểu hiện gen
# Tạo random seed để kết quả có thể tái tạo được
set.seed(123)

# Tạo dataframe với 8 gene và 8 mẫu (4 wt và 4 ko)
gene_data <- data.frame(
  wt1 = round(rnorm(8, mean = 5, sd = 1), 2),
  wt2 = round(rnorm(8, mean = 5, sd = 1), 2),
  wt3 = round(rnorm(8, mean = 5, sd = 1), 2),
  wt4 = round(rnorm(8, mean = 5, sd = 1), 2),
  ko1 = round(rnorm(8, mean = 3, sd = 1), 2),
  ko2 = round(rnorm(8, mean = 3, sd = 1), 2),
  ko3 = round(rnorm(8, mean = 3, sd = 1), 2),
  ko4 = round(rnorm(8, mean = 3, sd = 1), 2)
)

# Đặt tên cho các gene
rownames(gene_data) <- paste0("Gene", 1:8)

# In ra dataset
print(gene_data)

##        wt1  wt2  wt3  wt4  ko1  ko2  ko3  ko4
## Gene1 4.44 4.31 5.50 4.37 3.90 2.31 3.78 1.45
## Gene2 4.77 4.55 3.03 3.31 3.88 2.79 2.92 3.58
## Gene3 6.56 6.22 5.70 5.84 3.82 1.73 3.25 3.12
## Gene4 5.07 5.36 4.53 5.15 3.69 5.17 2.97 3.22
## Gene5 5.13 5.40 3.93 3.86 3.55 4.21 2.96 3.38
## Gene6 6.72 5.11 4.78 6.25 2.94 1.88 4.37 2.50
## Gene7 5.46 4.44 3.97 5.43 2.69 2.60 2.77 2.67
## Gene8 3.73 6.79 4.27 4.70 2.62 2.53 4.52 1.98

Thực hiện PCA

Hàm sẽ xem một cột là một chiều, sau khi PCA, thông tin về gen sẽ mất đi
Hàm trả về 3 giá trị: stdev, rotation và x. Trong đó x là tọa độ dữ liệu trong không gian mới và được sử dụng để vẽ

pca_result <- prcomp(t(gene_data), 
                     center = TRUE,  # Chuẩn hóa trung bình về 0
                     scale = TRUE)   # Chuẩn hóa độ lệch chuẩn về 1

Vẽ barplot

pca.var <- pca_result$sdev^2
pca.var.per <- pca.var/sum(pca.var)*100#phan tram của mỗi PCA var
barplot(pca.var.per)

Như vậy có thể thấy PCA 1 chủ yếu phân chia các nhóm

Tạo dataframe

PCA1 <- pca_result$x[,1]
PCA2 <- pca_result$x[,2]
sample <- row.names(pca_result$x)
group <-  substr(rownames(pca_result$x), 1, 2)
pca_data <- data.frame(sample, PCA1, PCA2, group)
  # Lấy 2 ký tự đầu (wt/ko)

Vẽ plot PCA

library(ggplot2)

## Warning: package 'ggplot2' was built under R version 4.2.3

library(dplyr)

## Warning: package 'dplyr' was built under R version 4.2.3

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

pca_data %>% ggplot(aes(x = PCA1, y = PCA2, color = group))+
    geom_point()+
    theme_bw()

PCA

huy

2024-10-26