diff --git a/R/calculate_lr.R b/R/calculate_lr.R
index 43d1da1f8fff6bad35723c15614c2bffce3ace90..72007f85aeafc7132f9f175f1e67565d4b343e9f 100644
--- a/R/calculate_lr.R
+++ b/R/calculate_lr.R
@@ -1,52 +1,81 @@
-#' Calculate lr values from summary file
+#' Calculate lr values from summary file and standardize lr values
 #'
-#' @param summary_file imported summary file (usually named 'summary')
+#' @param axiom_data imported Axiom files (obtained by calling read_axiom())
+#' @param standard whether or not to standardize the lr values ('none', 'mean' - the default , 'ref' or 'both')
+#' @param ref name of the sample to use for standardization (must be defined if method is "ref" or "both", defaults to NULL)
 #'
 #' @import dplyr
-#' @import tidyr
 #'
 #' @return A [tibble()].
 #' @export
 #'
 #' @examples
 #' \dontrun{
-#' calculate_lr(summary_file)
+#' calculate_lr(axiom_data, standard = "mean", ref = NULL)
 #' }
 
 calculate_lr <- function(summary_file) {
 
   # Set NULL variables
 
-  a <- b <- pivot_longer <- probeset_id <- NULL
-
-  # Extract A signal from summary file
-  message("Extracting A signal from summary file")
-  signal_a <- summary_file |>
-    dplyr::filter(!grepl(pattern = "NP", x = probeset_id)) |>
-    dplyr::filter(grepl("AX-[0-9]*-A", probeset_id)) |>
-    dplyr::mutate(probeset_id = stringr::str_remove(string = probeset_id, pattern = "-A")) |>
-    tidyr::pivot_longer(cols = -probeset_id,
-                        names_to = "file_name",
-                        values_to = "signal_a")
-
-  # Extract B signal from summary file
-  message("Extracting B signal from summary file")
-  signal_b <- summary_file |>
-    dplyr::filter(!grepl(pattern = "NP", x = probeset_id)) |>
-    dplyr::filter(grepl("AX-[0-9]*-B", probeset_id)) |>
-    dplyr::mutate(probeset_id = stringr::str_remove(string = probeset_id, pattern = "-B")) |>
-    tidyr::pivot_longer(cols = -probeset_id,
-                        names_to = "file_name",
-                        values_to = "signal_b")
-
-  # Calculate lr values
-  message("Calculating lr values")
-  lr <- signal_a |>
-    dplyr::left_join(signal_b) |>
-    dplyr::mutate(lr = log2(sqrt(signal_a^2 + signal_b^2)),
-                  .keep = "unused")
-
-  # Assign to Global Environment
-  assign(x = "lr", value = lr, pos = ".GlobalEnv")
+  # a <- b <- pivot_longer <- probeset_id <- NULL
+
+  # If standard = "none", calculate lr values
+  if (standard == "none") {
+    message("Calculating lr values")
+    lr <- axiom_data |>
+      dplyr::mutate(lr = log2(sqrt(signal_a^2 + signal_b^2)),
+                    .after = signal_b)
+  }
+
+  # If standard = "ref" or "both", check if ref is defined
+  if(standard  %in% c("ref", "both")){
+    print("Checking if reference sample is defined")
+    if(is.null(ref)){
+      stop("Reference sample is not defined")
+    }
+    if(!is.null(ref) & !ref %in% unique(axiom_data$file_name)) {
+      stop("Reference sample doesn't exist in the list of samples")
+    }
+  }
+
+  # If standard = "mean", calculate lr + lr_st_mean
+  if (standard == "mean") {
+    message("Calculating lr values and lr_st_mean values")
+    lr <- axiom_data |>
+      dplyr::mutate(lr = log2(sqrt(signal_a^2 + signal_b^2)),
+                    .after = signal_b) |>
+      dplyr::mutate(lr_mean = mean(lr, na.rm = TRUE),
+                    .by = probeset_id) |>
+      dplyr::mutate(lr_st_mean = lr - lr_mean,
+                    .after = lr) |>
+      dplyr::select(-lr_mean)
+  }
+
+  # If standard = "ref", calculatelr +  lr_st_ref
+  if (standard == "ref") {
+    message("Calculating lr values and lr_st_ref values")
+    lr <- axiom_data |>
+      dplyr::mutate(lr = log2(sqrt(signal_a^2 + signal_b^2)),
+                    .after = signal_b) |>
+      dplyr::mutate(lr_st_ref = lr - lr[file_name == ref],
+                    .by = probeset_id,
+                    .after = lr)
+  }
+
+  # If standard = "both", calculate lr + lr_st_ref & lr_st_mean
+  if (standard == "both") {
+    message("Calculating lr values, lr_st_mean and lr_st_ref values")
+    lr <- axiom_data |>
+      dplyr::mutate(lr = log2(sqrt(signal_a^2 + signal_b^2)),
+                    .after = signal_b) |>
+      dplyr::mutate(lr_st_ref = lr - lr[file_name == ref],
+                    lr_mean = mean(lr, na.rm = TRUE),
+                    .by = probeset_id) |>
+      dplyr::mutate(lr_st_mean = lr - lr_mean) |>
+      dplyr::select(file_name:lr, lr_st_mean, lr_st_ref, genotyping_call)
+  }
+
+  return(lr)
 
 }
diff --git a/R/standardise_lr.R b/R/standardise_lr.R
index f7b3495a1c94f87bab49229b6c37582c18b10d99..040110fd26d0e9259462544c7c51d0e8deeaab84 100644
--- a/R/standardise_lr.R
+++ b/R/standardise_lr.R
@@ -56,7 +56,4 @@ standardise_lr <- function(lr, method = "mean", ref = NULL) {
       dplyr::select(-lr_mean)
   }
 
-  # Assign to Global Environment
-  assign(x = "lr", value = lr, pos = ".GlobalEnv")
-
 }
diff --git a/man/calculate_lr.Rd b/man/calculate_lr.Rd
index d9a21914ecd766e173db37f3ed31c666d6827bbf..a9cf95d9d9785a3415b7938118e64c64aecd4cee 100644
--- a/man/calculate_lr.Rd
+++ b/man/calculate_lr.Rd
@@ -2,21 +2,25 @@
 % Please edit documentation in R/calculate_lr.R
 \name{calculate_lr}
 \alias{calculate_lr}
-\title{Calculate lr values from summary file}
+\title{Calculate lr values from summary file and standardize lr values}
 \usage{
 calculate_lr(summary_file)
 }
 \arguments{
-\item{summary_file}{imported summary file (usually named 'summary')}
+\item{axiom_data}{imported Axiom files (obtained by calling read_axiom())}
+
+\item{standard}{whether or not to standardize the lr values ('none', 'mean' - the default , 'ref' or 'both')}
+
+\item{ref}{name of the sample to use for standardization (must be defined if method is "ref" or "both", defaults to NULL)}
 }
 \value{
 A \code{\link[=tibble]{tibble()}}.
 }
 \description{
-Calculate lr values from summary file
+Calculate lr values from summary file and standardize lr values
 }
 \examples{
 \dontrun{
-calculate_lr(summary_file)
+calculate_lr(axiom_data, standard = "mean", ref = NULL)
 }
 }