limits of geom_histogram in a loop

[ZahraEM77]

Hi DS team- I hope you're well!

I am running a for loop which filters and plots my data according to 3 conditions (method, developmental stage, gene). I was wondering if there is a way to keep the scale of the histograms (the x and ylim) the same for all the plots of a single gene. E.g. all plots describing the expression level of gene x, (regardless of the developmental stage of the cell, or the method used), have the same limits. I'm hoping this will help me track how gene expression changes throughout different developmental stages.

This is my current code:

all_genes = unique(exonic_dataset_longer$Gene)
all_sorted = unique(exonic_dataset_longer$Sorted_Classification)
all_method = unique(exonic_dataset_longer$Rep)
library(ggplot2)
for (m in all_method ){
  exonic_data_m = filter(exonic_dataset_longer, Rep == m)
  for (s in all_sorted){
    exonic_data_m_s = filter(exonic_data_m, Sorted_Classification == s)
    
    for (g in all_genes){
      
      exonic_data_m_s_g = filter(exonic_data_m_s, Gene == g)
      p = ggplot(exonic_data_m_s_g, aes(Expr)) +
        geom_histogram(fill="royalblue", alpha=0.6) +
        ggtitle(paste0( g, " ", s, " ", m, ".png"))+
        theme(axis.text.x=element_text(angle=45, hjust=1, vjust=1))      
             plot = p, width = 6, height=6)
      
    }
    
  }}

Thank you so much for your insight and support- truly appreciate it!

[almenal]

Hi Zahra!

Sure, there is a way! The way you specify an x-axis range in ggplot is by using the family of functions scale_x_*, in this case you'll need to use scale_x_continuous(). Since you'll want the same limits for each gene, you can move the inner-most loop to the top, so that you can calculate the range of a gene first, and then apply it to all plots

In case it wasn't clear, your current loop is like this:

for (m in all_method ){
  ...
  for (s in all_sorted){
    ...
    for (g in all_genes){
      }
   }
}

so you should move the for (g in all_genes){ to the top. You can do the same thing you're doing of filtering dataframes sequentially to get a dataframe with only gene g.
Once you have that dataframe, the range of expression values will be from zero to the maximum (which you could keep in a variable called e.g. max_expr)
Then you can use those values to in scale_x_continuous(limits = c(0, max_expr)

for (g in all_genes){
  exonic_data_g = filter(exonic_data, Gene == g)
  max_expr = exonic_data_g$Expr

  for (s in all_sorted){
    ...
    for (m in all_method){

      exonic_data_m_s_g = ... 

      p = ggplot(exonic_data_m_s_g, aes(Expr)) +
        geom_histogram(fill="royalblue", alpha=0.6) +

        scale_x_continuous(limits = c(0, max_expr)) + # to keep the x axis constant for a given gene

        ggtitle(paste0( g, " ", s, " ", m, ".png"))+
        theme(axis.text.x=element_text(angle=45, hjust=1, vjust=1))      
             plot = p, width = 6, height=6)

      }
   }
}

Hope that solves your question!

[ZahraEM77]

Hi Almenal- Thank you so much for your insight- I really appreciate it!

I have implemented the changes to the code you suggested, and the x-axis limits of the various plots of a given gene have become consistent! However, now some of the genes don't show any expression, or the histogram is inaccurate relative to the data. The code I have used is the following:

all_genes = unique(exonic_dataset_longer$Gene)
all_sorted = unique(exonic_dataset_longer$Sorted_Classification)
all_method = unique(exonic_dataset_longer$Rep)
for (g in all_genes[1:4]){
exonic_data_g = filter(exonic_dataset_longer, Gene == g)
max_expr= exonic_data_g$Expr
for (s in all_sorted){
exonic_data_s_g = filter(exonic_data_g, Sorted_Classification == s)
for (m in all_method ){
exonic_data_m_s_g = filter(exonic_data_s_g, Rep == m)

  p = ggplot(exonic_data_m_s_g, aes(Expr)) +
    geom_histogram(fill="royalblue", alpha=0.6) +
    #scale_x_continuous(limits = c(0,max_expr))+
    ggtitle(paste0( g, " ", s, " ", m, ".png"))+
    theme(axis.text.x=element_text(angle=45, hjust=1, vjust=1))      
  ggsave(paste0("/Users/ZahraMusavi/Dropbox/ICL- MRes/Project 2/Veloctiy files/Exonic Selection Genes_pre-sel_sel-int/test/",
                g, " ", s, " ", m, ".png"),
         plot = p, width = 6, height=6)
  
}

}}

Could you please guide me as to what may be causing this?

Thank you again for your time!

[almenal]

Whoops, I forgot to write a max() over there. It should be

max_expr= max(exonic_data_g$Expr, na.rm = TRUE)

[ZahraEM77]

Thank you so so much- this is going to save me so much time! The graphs look great :)