Expression Quantitative Trait Loci (eQTL) Mapping

Module 6, Genetics & Genomics Winter School 2026

Author

Dr Solal Chauquet & Dr Zoe Hunter

Introduction

As we heard in Lecture 1, expression Quantitative Trait Loci (eQTL, or eQTLs) are genetic loci where variants, or single nucleotide polymorphisms (SNPs) at these loci are associated with increases or decreases in the expression of a specific gene. For example, figure 1 highlights how a SNP (in red) can be associated with gene expression. Here an A allele increases gene expression in a dose-response manner, with A homozygotes (AA, green curve) displaying higher levels of expression of a gene compared to G homozygotes (GG, black curve). Additionally, alleles are not expressed uniformly, with some variants rarer than others, as represented by the relative frequency of different genotypes on the y-axis.

Figure 1: Representation of an eQTL effect on a gene. Figure taken from Nica, A. C. & Dermitzakis, E. T. Expression quantitative trait loci: present and future. Philos. Trans. R. Soc. London. Ser. B, Biol. Sci. 368, 20120362 (2013).

Most eQTLs are found outside of coding regions and can be divided into two categories:

Cis-eQTls: found on the same chromosome as the gene they influence (usually within a 100,000 base pair, or 1Mb, window of that gene). Cis-eQTLs are thought to act on the gene directly through the regulation of enhancers, silencers, promoter regions or other regulatory elements of the gene.
Trans-eQTLs: found anywhere in the genome, further away from the gene they influence (> 500,000 base pairs, or 5Mb, away) or even on other chromosomes. They are thought to influence gene expression through regulation of biological pathways.

To investigate how to identify eQTLs, we shall split this practical into two parts:

Part 1 will start with simulating both genotype and expression data. This simulation approach will allow us to understand the structure of the data used for eQTL mapping, as well as investigate information relevant to QTL mapping more broadly.
Part 2 of the practical will explore the GTEx website and see how eQTLs can be used to investigate genome-wide association study (GWAS) results in real life.

Part 1: eQTL simulation

Connection to the cluster

To simulate the data for this practical, we will need to connect to the HPC cluster set for this class. You can reference the introduction to computing document to find how to connect to the HPC.

Alternatively, you can run the analysis on your local machine. Keep in mind that the simulation might take a while.

To begin with, we will set up the HPC folders to keep your analyses organised. The following bash script allows you to create the required three folders:

cd ~
mkdir eQTLPrac    # eQTL prac main folder
cd eQTLPrac       # Change directory to the eQTLPrac folder

# Inside the eQTLPrac folder, create the following subfolders to store results
mkdir Genotype    
mkdir Expression  
mkdir eQTL

Files created during the practical will be stored in these folders.

Gene expression & genotype data for QTL mapping

During the lecture, we saw the QTL mapping can be performed using a simple linear regression of the following form:

\[y = \beta_{0} + \beta_{1}x + \epsilon\]

With:

\(y\): Gene expression for each individual measured
\(\beta_{0}\): Intercept (mean gene expression across alleles)
\(\beta_{1}\): Effect of an allele on the gene expression
\(x\): Genotype value of each individual measured (0, 1, or 2)
\(\epsilon\): Error term of the model, or the part of \(y\) that is not explained by \(x\) (i.e. environment, batch effect)

To perform QLT mapping, we will have to simulate the gene expression, \(y\), as well as the genotype, \(x\), that will be used to infer the other parameters (\(\beta_{0}, \beta_{1}\)).

Genotype data

Genotypes are usually represented based on the number of alleles an individual carries at a specific locus. However, the allele used as a reference is arbitrary.

Taking this into account, we can represent individuals based on the number of alleles they carry. For example, the genotype of two individuals (C/G and G/G) at a specific locus can be shown as [1 0] or [1 2], depending on whether the C or G is used as the reference allele.

Question 1:

Think about a single genetic locus where the allele can be either A or T. How would you represent four individuals whose genotype at this locus are respectively:

Note: there are two possible sets of answers.

For the question above, we can represent the individuals based on the number of A alleles they carry: [2 1 1 0], or based on the number of T alleles they carry: [0 1 1 2].

Notice that we observe no difference in the annotation for homozygous individuals as they have the same number of A and T alleles. Similarly, the strand carrying the allele does not matter (for this kind of analysis), and A/T and T/A are both represented as 1.

To simulate genotypes, we therefore have to create several vectors containing 0, 1, or 2, representing the genotype of a single individual at different loci. Those vectors are then arranged into a matrix of the following form:

\[ \begin{bmatrix} 2 & 1 & 1 & 0 \\ 1 & 0 & 0 & 0 \\ 1 & 2 & 1 & 2 \\ 1 & 2 & 1 & 0 \end{bmatrix} \]

With the columns of the matrix representing each individual and the rows a specific genetic locus.

Simulation of genotypes

The following R code will now perform the simulation of genotypes. Start your R session and copy and paste the following into R:

# Load the libraries needed for the practical
library(tidyverse)
library(MASS)      
library(cowplot)
library(MetBrewer)

# Set seed makes sure we all get the same results (reproducibility)
set.seed(6543456)  

# Set the allele frequency and generate a vector of 0s, 1s, or 2s
frequency <- 0.5 # allele frequency
SNP <- rbinom(500, size = 2, frequency)

# Set the parameters for the genotype simulation
SNP_number  <- 1000
indv_number <- 500
p <- runif(SNP_number, min = 0, max = 1)

# Simulate the genotype matrix 
genotypes <- replicate(indv_number, rbinom(SNP_number, 2, p))
rownames(genotypes) <- paste0('SNP',  seq(1, nrow(genotypes)))
colnames(genotypes) <- paste0('Indv', seq(1, ncol(genotypes)))

# Inspect the resulting matrix
print(nrow(genotypes))
print(ncol(genotypes))
print(genotypes[1:10, 1:10])

Note: the set.seed function allows the code to be reproducible by fixing random processes (such as random sampling for simulation). A different seed would change the results of the downstream analysis.

Question 2:

How many SNPs were simulated?
How many individuals were simulated?
Given that the SNP3 reference allele is G and the alternative allele is C, what is the genotype of individual 5?

Exploration of allele frequency

Now that we have simulated the genotype data, we can calculate the frequency of the alleles simulated with the following code:

# Find the allele frequency of each SNP: average copies per individual, divided by 2
maf <- rowMeans(genotypes)/2 

# Then calculate the minor allele frequency from the mean frequency
maf <- pmin(maf, 1-maf)

# Create and save a histogram plot to visualise the distribution of the allele frequencies
jpeg('~/eQTLPrac/Genotype/HistogramMAFsimulated.jpeg',width = 21, height = 12, res = 300, units = 'cm', type = 'cairo')
truehist(maf, main = "Histogram of minor allele frequency", col = "light grey", nbins=100)
lines(density(maf), lty = 2, col = "dark red", lwd = 3)
dev.off()

If you are using the cluster to perform these analyses, you can download the plot that you just created by using the following command on your local machine (opening a new terminal window for downloading is always a great idea!):

scp <username>@203.101.xxx.xxx:~/eQTLPrac/Genotype/HistogramMAFsimulated.jpeg .

Question 3:

Look at the allele frequency of the genotype data you simulated:

What allele frequency is occurring the most?
In your opinion, why is allele frequency important for eQTL analysis? (Ask a tutor if your are not sure!)
What does the x-axis represent? Why is it limited to 0.5?

Allele frequency, particularly the minor allele frequency (MAF), is an important parameter during eQTL analysis due to the possible lack of representation of some genotypes. This is often the case for rare alleles where genotypes will often not be represented, affecting the statistical power of the analysis and/or resulting in false positives. For QTL analyses to be possible, we need to initially include individuals with all three genotype groups or to have at least two genotype groups present (i.e. [0 1], [1 2], or [0 2]).

Question 4:

Fill the table below using the genotype frequency derived from the Hardy-Weinberg principle for allele A, with a frequency \(p\) = 0.99.

What is the property of a linear regression that allows us to perform eQTL mapping when only 2 genotype groups are present?
Out of the three different populations in the table, which one(s) could be used to perform eQTL mapping for alleles with a frequency of 1%?

Table 1: Proportion of the possible genotypes for a genetic locus with a minor allele frequency q of 0.01 (calculated under the Hardy–Weinberg principle).

Genotype	Population genotype frequency	1,000 Individuals	10,000 Individuals	100,000 Individuals
AA Frequency: \(p^2\)
AT Frequency: \(2pq\)
TT Frequency: \(q^2\)

To identify eQTLs with a MAF of 1%, we would need a population larger than 10,000 individuals where we have access to two genotypes. With 198 heterozygous individuals expected in the population, we can therefore perform eQTL mapping. Given that the eQTL effect is expected to be additive (the effect of two alleles is twice the effect of one allele), the eQTL mapping can be performed with only 2 genotype groups.

The largest eQTL study to date, the eQTLgen Consortium, contains 31,684 individuals. This highlights that eQTL mapping for rare variants (MAF < 1%) is currently not possible. Allele frequency always needs to be kept in mind when performing eQTL mapping; exclusion of rare variants is currently always performed.

Gene expression data

Now that we have simulated genotype data, we need to create matching gene expression data. While gene expression is not normally distributed (RNA-sequencing and read-based sequencing technology generate discrete data that usually follows a negative binomial or Poisson distribution), most analyses will start by normalising the data. Simulating gene expression data can be performed either at the discrete level or at the normalised level.

For simplicity, this practical will simulate normally distributed data.

# Set the number of genes/QTLs
genesTotal      <- 1000
geneswithQTL    <- 50                          
geneswithoutQTL <- genesTotal - geneswithQTL

# Set the number of SNPs associated with genes
SNPs            <- rownames(genotypes)
SNPswithQTL     <- sample(SNPs, size = geneswithQTL)
SNPswithoutQTL  <- SNPs[-which(SNPs %in% SNPswithQTL)]


# Simulate the genes that DO have an eQTL

# Bulk RNA-seq mixes several cell types together (e.g. whole blood = many
# white-blood-cell types). We simulate 3 cell types per gene, then add them
# up to get the "bulk" value.

expMatrixAssociated <- do.call(cbind, lapply(SNPswithQTL, function(i) {

  # One mean per genotype (0, 1, 2 copies) per cell type
  meanCT1 <- c(rnorm(1, mean = 3), rnorm(1, mean = 5), rnorm(1, mean = 7))
  meanCT2 <- c(rnorm(1, mean = 3), rnorm(1, mean = 3), rnorm(1, mean = 3))
  meanCT3 <- c(rnorm(1, mean = 4), rnorm(1, mean = 2), rnorm(1, mean = 0))

  # genotypes[i, ] holds 0, 1, 2 for each person. Adding 1 turns those into
  # positions 1/2/3, so we pick the matching mean for that person's genotype,
  # then draw their expression around it.
  g <- genotypes[i, ] + 1
  df <- data.frame(
    ct1 = rnorm(indv_number, meanCT1[g]),
    ct2 = rnorm(indv_number, meanCT2[g]),
    ct3 = rnorm(indv_number, meanCT3[g])
  )

  # Bulk = sum of the three cell types for each person
  df$bulk <- rowSums(df, na.rm = TRUE)
  return(df)
}))

# Name the columns: Gene1_ct1, Gene1_ct2, Gene1_ct3, Gene1_bulk, Gene2_ct1, ...
colnames(expMatrixAssociated) <- paste0(
  'Gene', rep(1:geneswithQTL, each = 4), '_', colnames(expMatrixAssociated)
)


# Simulate the genes that do NOT have an eQTL

expMatrixNotAssociated <- do.call(cbind, lapply(SNPswithoutQTL, function(i) {

  meanForGene <- rnorm(1, mean = 10)

  df <- data.frame(
    ct1 = rnorm(indv_number, meanForGene),
    ct2 = rnorm(indv_number, meanForGene),
    ct3 = rnorm(indv_number, meanForGene)
  )
  df$bulk <- rowSums(df, na.rm = TRUE)
  return(df)
}))

colnames(expMatrixNotAssociated) <- paste0(
  'Gene', rep(1:geneswithoutQTL, each = 4), '_', colnames(expMatrixNotAssociated)
)


# Separate into a bulk dataframe and a cell-type-specific dataframe

# Bulk df: take only the 'bulk' columns
expMatrix <- cbind(
  expMatrixAssociated    %>% dplyr::select(contains('bulk')),
  expMatrixNotAssociated %>% dplyr::select(contains('bulk'))
)
colnames(expMatrix) <- paste0('Gene', 1:ncol(expMatrix))

# Cell-type df: take everything EXCEPT the 'bulk' columns
expMatrixCT <- cbind(
  expMatrixAssociated    %>% dplyr::select(!contains('bulk')),
  expMatrixNotAssociated %>% dplyr::select(!contains('bulk'))
)
colnames(expMatrixCT) <- paste0('Gene', rep(1:genesTotal, each = 3), '_', c('ct1', 'ct2', 'ct3'))

Question 5:

Run the code above and answer the following questions:

How many genes were simulated?
How many of those genes were associated with SNPs?
How many cell types were simulated?
How was the bulk expression created?
What omics technology does the bulk data correspond to?

The following code will generate plots showing the bulk expression level of one gene against the genotypes of several SNPs. You can change the code (by changing the name of the SNPs and genes) to visually inspect the association between different SNPs and genes.

### Test for SNPs:
SNPassociationPlot <- function(expMatrix, SNPID, GeneID) {
 ggplot(data.frame(snp = genotypes[SNPID,], y = expMatrix[,GeneID]),
 aes(x = factor(snp), y = y)) +
 ggtitle(paste0('Association between ', GeneID, ' and ', SNPID)) +
 geom_boxplot(fill = 'dark red') + geom_point(col = 'dark grey') +
 xlab("Reference allele count") +
 theme_minimal() + theme(plot.title = element_text(hjust = 0.5))
}

p1 <- SNPassociationPlot(expMatrix, SNPID = 'SNP182', GeneID = 'Gene2')
p2 <- SNPassociationPlot(expMatrix, SNPID = 'SNP243', GeneID = 'Gene2')
p3 <- SNPassociationPlot(expMatrix, SNPID = 'SNP921', GeneID = 'Gene2')
p4 <- SNPassociationPlot(expMatrix, SNPID = 'SNP564', GeneID = 'Gene2')
p <- cowplot::plot_grid(p1,p2,p3,p4)
ggsave(p, filename = '~/eQTLPrac/Expression/GenotypeVsExpressionPlot.jpeg', width = 14, height=14, dpi =
300)

If you are using the cluster, download the plot created using the following command:

scp <username>@203.101.xxx.xxx:~/eQTLPrac/Expression/GenotypeVsExpressionPlot.jpeg

While identifying SNP and gene expression pairs visually is already time consuming, the human genome is composed of 3.2 billion base pairs and roughly 20,000 genes, rendering it impossible to do manually. We need to use the linear regression that we previously described to filter the results based on significance.

eQTL mapping

Now that we have simulated both gene expression and genotypes, we will use linear regression to identify significant eQTLs.

Note that the plot function below will not work until you change the SNP ID value of the function.

# eQTL mapping for one gene
GeneID <- 'Gene1'

# Test each SNP for association with this gene's expression.
# For each SNP we fit a linear regression
Association <- data.frame()
for (SNPID in rownames(genotypes)) {
  test <- summary(lm(expMatrix[, GeneID] ~ genotypes[SNPID, ]))
  test <- as.data.frame(test$coefficients)[2, ]   
  rownames(test) <- SNPID
  Association <- rbind(Association, test)         
}

# Specify column names and reorder to show the most significant SNPs (smallest p-values) first
colnames(Association) <- c("Estimate", "Std.Error", "t_value", "P")
Association %>% arrange(P) %>% head() %>% print()

# Plot the most significant SNP
# Replace 'ChangeSNP' with the most significant SNP from the table above
signifQTL <- SNPassociationPlot(expMatrix, SNPID = 'ChangeSNP', GeneID = 'Gene1')
ggsave(signifQTL,
       filename = '~/eQTLPrac/Expression/significantBulkQTL.jpeg',
       width = 7, height = 7, dpi = 300)

Question 7:

What SNP (if any) is significantly associated with Gene 1?
Fill the gap in the code (SNPID) with the significantly associated SNP and investigate the association visually.
Modify the previous code to investigate associations with gene 982.

The association test between the gene and 1,000 SNPs for 500 individuals that we just performed took only a few seconds. However, this toy example does not represent the scale of the human genome with its more than 3 billion base pairs. eQTL analyses quickly result in prohibitive computation time as we increase the number of SNPs and individuals tested.

You can modify the variables SNP_number and indv_number and rerun the code up to this point to simulate a larger number of individuals an SNPs. However, increasing the number of individuals to 2,000 will already heavily increase the computational time needed to run this (poorly optimised) script.

Software such as MatrixQTL and fastQTL have been developed to decrease the computational resources and time necessary for eQTL analyses, While we will not go into details on the inner working here, the underlying mechanisms of that software remain similar to the analysis performed within this practical (linear regression - the backbone of genetics). Methodology used to improve computational efficiency ranges from limiting the SNPs tested for a gene to the closest SNPs, to developing mathematical approximations to computationally heavy calculations.

eQTL mapping with interaction

Our simulation of gene expression data was based on the presence of three different cell types measured, which is representative of bulk-RNA sequencing. We will now see what those QLTs look like when we decompose them across cell type.

Question 8:

Modify and run the script below to visually inspect Gene 1 with the SNP that you previously found to be significantly associated.
Does the eQTL identified previously represent a specific cell type?

set.seed(58944)
# Using the cell-type specific expression matrix 
expMatrixCT <- expMatrixCT %>% mutate(Indv = paste0('Indv', seq(1, indv_number)))
expMatrixCTlonger <- expMatrixCT %>% pivot_longer(cols = -Indv,
                                                  names_to = 'geneID',
                                                  values_to = 'expression') %>%
  mutate(
    gene = str_split(geneID, '_', simplify = T)[, 1],
    cellType = str_split(geneID, '_', simplify = T)[, 2]
  )
genotypeToMerge <- t(genotypes) %>% as.data.frame() %>%
  mutate(Indv = rownames(.))
expMatrixCTlonger <- left_join(expMatrixCTlonger, genotypeToMerge, by = 'Indv')


# Plot gene 1:
# Change the code in Red to inspect the gene of interest:
ggplot(
  expMatrixCTlonger %>% filter(gene == 'Gene1'),
  aes(x = factor(SNPID), y = expression, fill = cellType)
) +
  geom_point(position = position_jitterdodge()) +
  facet_wrap( ~ cellType) +
  geom_boxplot() +
  theme_minimal() +
  scale_fill_manual(values = met.brewer(n = 3, 'Hokusai1')) -> interacionPlot
ggsave(
  interacionPlot,
  filename = '~/eQTLPrac/Expression/interacionQTL.jpeg',
  width = 7,
  height = 7,
  dpi = 300
)

As we can see, the eQTL observed previously was produced by a difference in gene expression between the different cell types. We do not observe any eQTL when the cell type information is known.

We will now investigate eQTLs when different cell types are included. Our previous linear regression can be written with interaction between cell type and genotype as follows:

\[y = \beta_{0} + \beta_{1}x_{1} + \beta_{2}x_{2} + \beta_{3}x_{1}x_{2} + \epsilon\]

With:

\(y\): Gene expression for each individual measured
\(\beta_{0}\): Intercept (mean gene expression across alleles)
\(\beta_{1}\): Effect of an allele on the gene expression
\(x\): Genotype value of each individual measured (0, 1, or 2)
\(\beta_{2}\): Overall effect of cell type on gene expression
\(x_{2}\): Cell type of origin
\(\beta_{3}\): Effect of genotype on gene expression within different cell types
\(\epsilon\): Error term of the model

Using this model, our term of interest will be \(\beta_{3}\), representing the effect of one cell type compared to the overall effect of the genotype.

Question 9:

Run the code below and modify the plotting function to output a significant interaction for Gene 999.

What is the interaction observed between cell type, your significant SNP, and Gene 999?
What does the bulk QTL data look like for your significant SNP and Gene 999? Is it a significant SNP?

interactionResults <- data.frame()
for (snp in SNPs) {
 test <- expMatrixCTlonger %>% filter(gene == 'Gene985') %>%
 dplyr::select(expression, gene, cellType, snp) %>% mutate(variant = snp)
 colnames(test) <- c('expression', 'gene', 'cellType', 'genotype', 'variant')

 lmTest <- broom::tidy(summary(lm(expression ~ genotype +
 cellType + genotype*cellType,
 data = test)))
 lmTest$SNP <- snp
 interactionResults <- rbind(interactionResults, lmTest)
}
interactionResults %>%
 filter(str_detect(term, ':')) %>%
 arrange(p.value) %>% head()
# Change the SNP value in Red to the most significant SNP:
ggplot(expMatrixCTlonger %>% filter(gene == 'Gene985'),
 aes(x = factor(SNP), y = expression, fill = cellType)) +
 geom_point(position = position_jitterdodge()) +
 facet_wrap(~cellType) +
 geom_boxplot() +
 theme_minimal() +
 scale_fill_manual(values = met.brewer( n = 3, 'Hokusai1')) -> signifInteracionQLT
ggsave(signifInteracionQLT, filename = '~/eQTLPrac/Expression/significantInteractionQTL.jpeg', width = 7, height=7, dpi = 300)

Part 2: Real world eQTL mapping

Genotype-Tissue Expression Quantitative Trait Loci

We will now investigate real-world eQTL data. If you do not have the time to go through the whole document, an answer sheet will be provided at the end of the practicals to go through in your own time.

For this part, we will go to the GTEx website. You can access it through this here

The GTEx consortium collected post-mortem samples for 948 donors. We know that eQTLs are dynamic and evolve over time and with exposure to the environment. Characteristics such as sex, age or disease status can influence eQTL associations and are therefore important.

Question 10:

On the GTEx website, look for the sample characteristics that could influence an eQTL association study.

Hint: Navigate to the Tissue & Sample statistics page.

eQTLs are influenced by both age, sex, and ancestries; the observed unbalanced number of male and females as well as a largely white and aging (84.6% white, 68.1% of samples older than 50) cohort, therefore, need to be taken into account when performing eQTL analysis. Additionally, the cohort can be split in half with younger donors succumbing to traumatic injury while older donors displaying non-traumatic pathologies.

Sample characteristics, therefore, need to be considered when performing QTL mapping. You can read the landmark GTEx publication to observe which sample characteristics were corrected for when testing for QTL associations.

Investigation of GWAS signal

We will now investigate a real example of an eQTL association. For this, we will start by looking at a genome-wide association study of lipids published in 2013 in nature genetics:

Discovery and refinement of loci associated with lipid levels

Question 11:

Read the abstract of the GWAS paper - what is the goal of this paper?

This paper aimed to identify the genetic control of blood lipid levels. As such, they identified associations between SNPs and blood lipid levels. They then mapped those SNPs to the closest genes, concluding that this means these genes had a role in blood lipid levels.

We will investigate how eQTL analysis can give us more information regarding the genetic control of blood low-density lipoprotein (LDL) cholesterol.

Open the supplementary figures from the paper and go to the supplementary table 3.

Question 12:

Find the gene with the strongest negative effect on LDL blood levels
What is the impact of each alternative allele (i.e. look at the effect column)?
If the average person has an LDL blood level of 209.7mg/dL, what would be the expected LDL level of an individual with a genotype of G/G, G/T, or T/T at locus rs6511720?

Let’s investigate the effect of rs6511720, the SNP associated with the highest decrease in LDL blood levels. Search the GTEx website for rs6511720 and answer the following questions.

Question 13:

With which genes is rs6511720 associated?
In which tissues are those associations located?
Do you think that a change in gene expression could be responsible for the association observed between LDL levels and rs6511720?

We will now look at genetic loci associated with LDL cholesterol levels. rs12916 is associated with HMGCR, a gene coding for HMG-CoA reductase, an enzyme playing a central role in cholesterol synthesis. Let’s investigate eQTLs associated with rs12916, and search the GTEx website for rs12916.

Question 14:

In which tissue is rs12916 associated with HMGCR?
Where does the SNP fall? (hint: open the IGV browser)
Do you think that a change in gene expression could be responsible for the association observed between LDL levels and rs12916?

In conclusion, eQTLs can help with interpreting the functional significance of GWAS signals. They can provide biological interpretation of non-coding variants, helping to hint at the mechanisms underlying complex traits and diseases.

Part 3: Extension

This part is completely optional.

Here are a few extension questions to help understand what QTLs are and how we can investigate them.

Extension question 1:

The code used to generate bulk RNA-seq at the start of this practical assumes an equal mixture of all three cell types.

Do you think that it is realistic?
Does the same mixture of cell types apply to all genes?

This underlying assumption is not biologically accurate. We made it here since the goal of this practical is to gain a better understanding of how QTL studies work, not to generate biologically accurate data.

Extension question 2:

Modify the code generating bulk gene expression data to include reads coming from different proportions of cell types.
How does this influence the identification of QTL at the bulk or single-cell level?

Good luck, this one is harder!