monitor_relion_classification.R

#! /usr/bin/env Rscript --vanilla

# For each class of a RELION Class2D or Class3D job, plot the number of particles
# as a function of the iteration number. This is useful to monitor convergence:
# do classes still significantly vary in size after a given number of iterations?
# Uses the summary tsv files generated by count_particles.awk.

# Load helpful packages
library(tidyverse)

# List summary files and load them in dataframes
my_files <- list.files(path = ".",
                       pattern = "_data.star.tsv")
my_datasets <- map(.x = my_files,
                   .f = read_tsv,
                   skip = 2)
names(my_datasets) <- my_files

# For each dataframe, add a column containing the iteration number
iteration_numbers <- seq_along(my_datasets) - 1
my_datasets <- map2(.x = my_datasets,
                    .y = iteration_numbers,
                    .f = ~ mutate(.data = .x, iteration = as.integer(.y)))

# Plot results based on particle number
my_datasets %>%
    bind_rows() %>%
    mutate(Class = as.factor(Class)) %>%
    ggplot(aes(x = iteration, y = Particles, color = Class)) +
    geom_line() +
    theme_bw() +
    xlab("Iteration") +
    ylab("Particles")

# Save plot
ggsave(filename = "relion_classification_progress_particles.png",
       device = "png",
       units = "cm",
       width = 24.27,
       height = 15,
       dpi = 300)

# Plot results based on class distribution
my_datasets %>%
    bind_rows() %>%
    mutate(Class = as.factor(Class)) %>%
    ggplot(aes(x = iteration, y = `% of total`, color = Class)) +
    geom_line() +
    theme_bw() +
    xlab("Iteration") +
    ylab("% of total")

# Save plot
ggsave(filename = "relion_classification_progress_distribution.png",
       device = "png",
       units = "cm",
       width = 24.27,
       height = 15,
       dpi = 300)

file.remove("Rplots.pdf")

# Clear environment
remove(list = ls())