Gene Expression in an Alzheimer’s Patient
Liam Garvey
2026-09-15
Purpose
The purpose of this project is to use exploratory data analysis to examine gene expression in brain tissue from a 90-year-old woman in the RUSH Alzheimer’s study. Three plots will be used to look at the distribution of gene expression, identify the most highly expressed genes, and examine the relationship between gene length and expression.
# Load tidyverse for organizing data and creating plots.
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.2.1 ✔ readr 2.2.0
## ✔ forcats 1.0.1 ✔ stringr 1.6.0
## ✔ ggplot2 4.0.3 ✔ tibble 3.3.1
## ✔ lubridate 1.9.5 ✔ tidyr 1.3.2
## ✔ purrr 1.2.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors# Import ENCODE data ------------------------------------
# Read the file into R and save as EncodeData
EncodeData <- read.table("ENCFF166SFX.tsv", header = TRUE)
# Show first 6 rows, display column names, show structure, and give a summary of the data
head(EncodeData)## gene_id transcript_id.s. length effective_length expected_count TPM FPKM
## 1 10904 10904 93 0 0 0 0
## 2 12954 12954 94 0 0 0 0
## 3 12956 12956 72 0 0 0 0
## 4 12958 12958 82 0 0 0 0
## 5 12960 12960 73 0 0 0 0
## 6 12962 12962 72 0 0 0 0
## posterior_mean_count posterior_standard_deviation_of_count pme_TPM pme_FPKM
## 1 0 0 0 0
## 2 0 0 0 0
## 3 0 0 0 0
## 4 0 0 0 0
## 5 0 0 0 0
## 6 0 0 0 0
## TPM_ci_lower_bound TPM_ci_upper_bound TPM_coefficient_of_quartile_variation
## 1 0 0 0
## 2 0 0 0
## 3 0 0 0
## 4 0 0 0
## 5 0 0 0
## 6 0 0 0
## FPKM_ci_lower_bound FPKM_ci_upper_bound
## 1 0 0
## 2 0 0
## 3 0 0
## 4 0 0
## 5 0 0
## 6 0 0
## FPKM_coefficient_of_quartile_variation
## 1 0
## 2 0
## 3 0
## 4 0
## 5 0
## 6 0names(EncodeData)## [1] "gene_id"
## [2] "transcript_id.s."
## [3] "length"
## [4] "effective_length"
## [5] "expected_count"
## [6] "TPM"
## [7] "FPKM"
## [8] "posterior_mean_count"
## [9] "posterior_standard_deviation_of_count"
## [10] "pme_TPM"
## [11] "pme_FPKM"
## [12] "TPM_ci_lower_bound"
## [13] "TPM_ci_upper_bound"
## [14] "TPM_coefficient_of_quartile_variation"
## [15] "FPKM_ci_lower_bound"
## [16] "FPKM_ci_upper_bound"
## [17] "FPKM_coefficient_of_quartile_variation"str(EncodeData)## 'data.frame': 59526 obs. of 17 variables:
## $ gene_id : chr "10904" "12954" "12956" "12958" ...
## $ transcript_id.s. : chr "10904" "12954" "12956" "12958" ...
## $ length : num 93 94 72 82 73 72 74 82 73 73 ...
## $ effective_length : num 0 0 0 0 0 0 0 0 0 0 ...
## $ expected_count : num 0 0 0 0 0 0 0 0 0 0 ...
## $ TPM : num 0 0 0 0 0 0 0 0 0 0 ...
## $ FPKM : num 0 0 0 0 0 0 0 0 0 0 ...
## $ posterior_mean_count : num 0 0 0 0 0 0 0 0 0 0 ...
## $ posterior_standard_deviation_of_count : num 0 0 0 0 0 0 0 0 0 0 ...
## $ pme_TPM : num 0 0 0 0 0 0 0 0 0 0 ...
## $ pme_FPKM : num 0 0 0 0 0 0 0 0 0 0 ...
## $ TPM_ci_lower_bound : num 0 0 0 0 0 0 0 0 0 0 ...
## $ TPM_ci_upper_bound : num 0 0 0 0 0 0 0 0 0 0 ...
## $ TPM_coefficient_of_quartile_variation : num 0 0 0 0 0 0 0 0 0 0 ...
## $ FPKM_ci_lower_bound : num 0 0 0 0 0 0 0 0 0 0 ...
## $ FPKM_ci_upper_bound : num 0 0 0 0 0 0 0 0 0 0 ...
## $ FPKM_coefficient_of_quartile_variation: num 0 0 0 0 0 0 0 0 0 0 ...summary(EncodeData)## gene_id transcript_id.s. length effective_length
## Length :59526 Length :59526 Min. : 8 Min. : 0
## N.unique :59526 N.unique :59526 1st Qu.: 387 1st Qu.: 287
## N.blank : 0 N.blank : 0 Median : 808 Median : 708
## Min.nchar: 4 Min.nchar: 4 Mean : 1339 Mean : 1240
## Max.nchar: 24 Max.nchar: 3455 3rd Qu.: 1823 3rd Qu.: 1723
## Max. :205012 Max. :204912
## expected_count TPM FPKM
## Min. : 0.0 Min. : 0.00 Min. : 0.00
## 1st Qu.: 0.0 1st Qu.: 0.00 1st Qu.: 0.00
## Median : 0.0 Median : 0.00 Median : 0.00
## Mean : 339.3 Mean : 16.80 Mean : 19.55
## 3rd Qu.: 46.0 3rd Qu.: 1.52 3rd Qu.: 1.77
## Max. :1251092.4 Max. :267430.21 Max. :311258.58
## posterior_mean_count posterior_standard_deviation_of_count pme_TPM
## Min. : 0.0 Min. : 0.000 Min. : 0.00
## 1st Qu.: 0.0 1st Qu.: 0.000 1st Qu.: 0.09
## Median : 0.0 Median : 0.000 Median : 0.30
## Mean : 339.3 Mean : 1.965 Mean : 16.80
## 3rd Qu.: 46.0 3rd Qu.: 0.000 3rd Qu.: 3.00
## Max. :1068373.8 Max. :22576.980 Max. :217493.25
## pme_FPKM TPM_ci_lower_bound TPM_ci_upper_bound
## Min. : 0.00 Min. :0.000e+00 Min. :0.000e+00
## 1st Qu.: 0.11 1st Qu.:0.000e+00 1st Qu.:2.530e-01
## Median : 0.36 Median :3.000e-03 Median :6.790e-01
## Mean : 20.22 Mean :1.502e+01 Mean :1.923e+01
## 3rd Qu.: 3.61 3rd Qu.:1.072e+00 3rd Qu.:4.465e+00
## Max. :261833.14 Max. :2.169e+05 Max. :2.181e+05
## TPM_coefficient_of_quartile_variation FPKM_ci_lower_bound FPKM_ci_upper_bound
## Min. :0.0000 Min. :0.000e+00 Min. :0.000e+00
## 1st Qu.:0.0759 1st Qu.:0.000e+00 1st Qu.:3.040e-01
## Median :0.3873 Median :4.000e-03 Median :8.180e-01
## Mean :0.3677 Mean :1.808e+01 Mean :2.316e+01
## 3rd Qu.:0.6553 3rd Qu.:1.292e+00 3rd Qu.:5.370e+00
## Max. :0.7098 Max. :2.612e+05 Max. :2.625e+05
## FPKM_coefficient_of_quartile_variation
## Min. :0.00000
## 1st Qu.:0.07591
## Median :0.38733
## Mean :0.36768
## 3rd Qu.:0.65535
## Max. :0.70981The ENCODE data is stored as a data frame containing gene identifiers and measurements of gene expression. Each row is one gene, while the columns are gene length, expected count, TPM, and FPKM.
Figure 1: Gene-Length Distribution
The first figure shows gene length as gene size can affect the number of sequencing fragments assigned to a gene. Longer genes have more bases where sequencing fragments can align which can skew the analysis.
# Prepare gene-length data ------------------------------
# Keep genes longer than zero but shorter than 10,000 base pairs.
ShortGenes <- EncodeData %>%
filter(
length > 0,
length < 10000
)
# Plot gene-length distribution -------------------------
# Create a histogram showing the distribution of gene lengths.
GeneLength <- ggplot(
ShortGenes,
aes(x = length)
) +
# Group gene lengths into 250-base-pair intervals.
geom_histogram(
binwidth = 250,
fill = "green",
color = "white"
) +
# Add a title and axis labels.
labs(
title = "Distribution of Gene Lengths",
subtitle = "Genes shorter than 10,000 base pairs",
x = "Gene Length in Base Pairs",
y = "Number of Genes"
)
# Display the completed graph.
GeneLength
Figure 1 shows the distribution of genes that are shorter than 10,000 base pairs. Most of the genes included in the graph are found toward the shorter end of this range, while fewer genes are closer to 10,000 base pairs. Gene length is important when looking at expression data as longer genes have more locations for sequencing reads to line up which can cause longer genes to have higher expected counts.
Figure 2: Expected Count Compared with FPKM
Expected count estimates the number of sequencing fragments assigned to each gene. FPKM normalizes expression based on both gene length and the total number of sequenced fragments.
# Prepare expression data -------------------------------
# Keep genes with positive expected count and FPKM values.
ExpressedGenes <- EncodeData %>%
filter(
expected_count > 0,
FPKM > 0
)
# Compare expected count with FPKM -----------------------
# Create a scatterplot comparing expected count and FPKM.
CountvFPKM <- ggplot(
ExpressedGenes,
aes(
x = expected_count,
y = FPKM,
color = length
)
) +
# Add one point for each expressed gene.
geom_point(
alpha = 0.4,
size = 1
) +
# Use logarithmic scales because both measurements contain a large ranges of values.
scale_x_log10() +
scale_y_log10() +
# Apply a color gradient based on gene length.
scale_color_gradient(
low = "gold",
high = "darkred"
) +
# Add a line showing the overall relationship.
geom_smooth(
method = "lm",
color = "black",
se = FALSE
) +
# Add title and axis labels.
labs(
title = "Expected Gene Count Compared with FPKM",
subtitle = "Each point represents one expressed gene",
x = "Expected Count on Log10 Scale",
y = "FPKM on Log10 Scale",
color = "Gene Length"
)
# Display the completed graph.
CountvFPKM## `geom_smooth()` using formula = 'y ~ x'
Figure 2 shows the relationship between expected count and FPKM. As expected count increases, FPKM generally increases as well. Tthe points do not form a perfect line as FPKM accounts for gene length and the total number of fragments in the samples. The colors show the lengths of the genes and help show why genes with similar expected counts can have different FPKM values.
Figure 3: Genes with the Highest FPKM
The final figure uses FPKM to identify the 15 genes with the highest normalized expression in this brain sample.
# Select genes with the highest FPKM --------------------
# Arrange genes from highest to lowest FPKM and keep the first 15 genes.
TopGenes <- EncodeData %>%
arrange(desc(FPKM)) %>%
slice_head(n = 15)
# Plot the highest FPKM values ---------------------------
# Create a bar graph of the most highly expressed genes.
TopFPKM <- ggplot(
TopGenes,
aes(
x = reorder(gene_id, FPKM),
y = FPKM,
fill = FPKM
)
) +
# Create one bar for each gene, turn the graph so the gene IDs are easier to read and use a log scale because the FPKM values have a large range.
geom_col() +
coord_flip() +
scale_y_log10() +
# Color the bars according to their FPKM values.
scale_fill_gradient(
low = "skyblue",
high = "navy"
) +
# Add title and axis labels.
labs(
title = "Genes with Highest FPKM Values",
x = "Ensembl Gene ID",
y = "FPKM on Log10 Scale",
fill = "FPKM"
)
# Display the completed graph.
BlueGenes <- TopFPKM
BlueGenes
Figure 3 shows that a small number of genes have much higher normalized expression than other genes. These genes produced some of the most abundant RNA transcripts in this brain sample. Looking up their Ensembl gene IDs could help determine whether they are involved in brain function, energy production, protein synthesis, or processes linked to Alzheimer’s disease.