Load required libraries
install.packages(c("MetaCycle", "circular"))
library(MetaCycle)
library(circular)
Read gene expression data from CSV file
expression_data <- read.csv("C:/My own analysis/R-trials/rees_time_trial.csv", row.names = 1)
Remove genes with 0 TPM at all time points
expression_data <- expression_data[rowSums(expression_data != 0, na.rm = TRUE) > 0, ]
Extract timepoints from column names
timepoints <- as.numeric(gsub("X", "", colnames(expression_data)))
Check for any NA values in timepoints
if (any(is.na(timepoints))) {
stop("The 'timepoints' contain NA or NaN values.")
}
Write data to temporary text file in the required format
temp_file <- tempfile(fileext = ".txt")
write.table(expression_data, temp_file, sep = "\t", quote = FALSE, col.names = NA)
Perform MetaCycle analysis with adjusted parameters
mc_result <- meta2d(temp_file, time = timepoints, filestyle = "txt", minper = 12, maxper = 35, adjustPhase = "predictedPer")
Extract cycling genes
cycling_genes <- mc_result$cycling.genes
Calculate circadian phase (CT) for each transcript
phase_estimates <- mc_result$phase.estimates
period_estimates <- mc_result$period.estimates
circadian_phase <- (phase_estimates * 24) / period_estimates
Calculate circular phase means
circular_phase_means <- circular::mean.circular(circular_phase, type = "angle", na.rm = TRUE)
Print cycling genes
print(cycling_genes)
Print circadian phase estimates
print(circadian_phase)
Print circular phase means
print(circular_phase_means)