SciCodePile/SciCode-Domain-Code
DATA1: Domain-Specific Code Dataset Dataset Overview DATA1 is a large-scale domain-specific code dataset focusing on code samples from interdisciplinary fields such as biology, chemistry, materials science, and related areas. The dataset is collected and organized from GitHub repositories, covering 178 different domain topics with over 1.1 billion lines of code. Dataset Statistics Total Datasets: 178 CSV files Total Data Size: ~115 GB Total Lines… See the full description on the dataset page: https://huggingface.co/datasets/SciCodePile/SciCode-Domain-Code.
42.4k
1"keyword","repo_name","file_path","file_extension","file_size","line_count","content","language"
2"Homeostasis","MathOnco/EvoFreq","img/readme.img.code.R",".R","3345","75","library(gridExtra)3 4# Note: this can be copy and pasted after installed5 6data(""example.easy.wide"") # Load a simple Data Frame example7str(example.easy.wide) # Inspect the data structure8 9# 'data.frame': 8 obs. of 10 variables:10# $ parents: num 0 1 1 3 1 5 5 511# $ clones : num 1 2 3 4 5 6 7 812# $ 1 : num 1 0 0 0 0 0 0 013# $ 2 : num 100 5 0 0 0 0 0 014# $ 3 : num 200 100 5 0 0 0 0 015# $ 4 : num 400 5 100 1 1 0 0 016# $ 5 : num 0 0 200 100 100 1 0 117# $ 6 : num 0 0 200 125 200 10 1 1518# $ 7 : num 0 0 300 200 300 20 10 2519# $ 8 : num 0 0 300 300 300 25 25 10020 21# You have A column of parents and a column of clones then you have a column for each of the timepoints with sizes for that clone.22 23# Then get the frequency data. (Use ?get_freq_dynamics for options)24clone_dynamics_df <- get_freq_dynamics(example.easy.wide[,seq(3,10)], example.easy.wide$clones, example.easy.wide$parents, clone_cmap = ""magma"")25 26# Create the plot (shown on the left below)27evo_freq_p <- plot_freq_dynamics(clone_dynamics_df)28print(evo_freq_p)29 30# We can also choose to update the colors or do this during the first creation. (shown on the right below)31clone_dynamics_df_jet <- update_colors(clone_dynamics_df, clone_cmap = ""jet"")32evo_freq_p_jet <- plot_freq_dynamics(clone_dynamics_df_jet)33print(evo_freq_p_jet)34 35p <- arrangeGrob(evo_freq_p, evo_freq_p_jet, ncol=2, nrow=1)36ggsave(plot=p, filename=""img/easy.wide.image.png"", units=c(""mm""),width=317.5, height=79.375, dpi=300, device=""png"")37 38# Note: this can be copy and pasted after installed39 40# In this example you have two files. One is the edge list of clones and their parents41# The other file is the sizes over time for those clones42# examine the structures to see how to format43# example.easy.long.edges now has an attribute column as well44data(""example.easy.long.edges"")45data(""example.easy.long.sizes"")46 47# Use the long_to_wide_size_df function to get the right data structure.48wide_df <- long_to_wide_size_df(long_pop_sizes_df = example.easy.long.sizes, time_col_name = ""Time"", clone_col_name = ""clone"", parent_col_name = ""parent"", size_col_name = ""Size"", edges_df = example.easy.long.edges)49clones <- wide_df$clones50parents <- wide_df$parents51size_df <- wide_df$wide_size_df52 53clone_dynamics_df <- get_freq_dynamics(size_df, clones, parents, clone_cmap = ""inferno"")54evo_freq_p <- plot_freq_dynamics(clone_dynamics_df)55 56# Add custom ggplot features57evo_freq_labeled_p <- evofreq_clone_labels(clone_dynamics_df, clone_list=c(3,5,6), extant_only=F, evo_freq_p=evo_freq_p, apply_labels=T)58 59p <- arrangeGrob(evo_freq_p, evo_freq_labeled_p, ncol=2, nrow=1)60ggsave(plot=p, filename=""img/easy.long.image.png"", units=c(""mm""),width=317.5, height=79.375, dpi=300, device=""png"")61 62# Network plots63data(""example.wide"")64time_col_idx <- which(substr(colnames(example.wide), start=1, stop = 1)==""X"")65attribute_col_idx <- which(substr(colnames(example.wide), start=1, stop = 1) != ""X"")66colnames(example.wide)[time_col_idx] <- substring(colnames(example.wide), first=2)[time_col_idx]67attribute_df <- example.wide[, attribute_col_idx]68size_df <- example.wide[, time_col_idx]69parents <- example.wide$parent70clones <- example.wide$clone71 72clone_network <- get_clone_network(size_df, clones, parents)73evofreq_plot <- plot_clone_network(clone_network, layout=""circlepack"")74 75 76","R"
77"Homeostasis","MathOnco/EvoFreq","R/example.easy.long.edges.R",".R","418","16","#' Clone edges in Easy Long Format78#'79#' A dataset for easy understanding of EvoFreq data structures.80#'81#' @format A data frame with 6 rows and 3 variables:82#' \describe{83#' \item{parent}{Parent of clone}84#' \item{clone}{Clone}85#' \item{fitness}{Clone attribute that to color by}86#' }87#' @source Custom Made.88#' @examples89#' ?example.easy.long.edges90#' data(""example.easy.long.edges"")91""example.easy.long.edges""92","R"
93"Homeostasis","MathOnco/EvoFreq","R/Example.Wide.R",".R","2073","56","#' Clone information dataset in Wide Format94#'95#' A dataset from a simulation taken from ""Niche engineering drives early passage through an immune bottleneck in progression to colorectal cancer""96#'97#' @format A data frame with 6591 rows and 43 variables:98#' \describe{99#' \item{master_random_seed}{Simulation metadata}100#' \item{macrophage_suppression}{Simulation metadata}101#' \item{macrophage_growth_benefit}{Simulation metadata}102#' \item{pdl1_protection}{Simulation metadata}103#' \item{tumor_random_seed}{Simulation metadata}104#' \item{tumor_lineage}{Simulation metadata}105#' \item{tumor_status}{Simulation metadata}106#' \item{onset}{Simulation metadata}107#' \item{parent}{Simulation metadata}108#' \item{og_parent}{Simulation metadata}109#' \item{clone}{Simulation metadata}110#' \item{type}{Simulation metadata}111#' \item{division_rate}{Simulation metadata}112#' \item{n_driver_mutations}{Simulation metadata}113#' \item{new_antigenicity}{Simulation metadata}114#' \item{og_antigenicity}{Simulation metadata}115#' \item{final_antigenicity}{Simulation metadata}116#' \item{macrophage}{Simulation metadata}117#' \item{pdl1}{Simulation metadata}118#' \item{immune_susceptibility}{Simulation metadata}119#' \item{detected}{Simulation metadata}120#' \item{size_at_time_of_detection}{Simulation metadata}121#' \item{X2227}{Timepoint}122#' \item{X2257}{Timepoint}123#' \item{X2287}{Timepoint}124#' \item{X2317}{Timepoint}125#' \item{X2347}{Timepoint}126#' \item{X2377}{Timepoint}127#' \item{X2407}{Timepoint}128#' \item{X2437}{Timepoint}129#' \item{X2467}{Timepoint}130#' \item{X2497}{Timepoint}131#' \item{X2527}{Timepoint}132#' \item{X2557}{Timepoint}133#' \item{X2587}{Timepoint}134#' \item{X2617}{Timepoint}135#' \item{X2647}{Timepoint}136#' \item{X2677}{Timepoint}137#' \item{X2707}{Timepoint}138#' \item{X2737}{Timepoint}139#' \item{X2767}{Timepoint}140#' \item{X2797}{Timepoint}141#' \item{X2808}{Timepoint}142#' }143#' @source \url{https://www.biorxiv.org/content/10.1101/623959v2}144#' @examples145#' if (require(""EvoFreq"")) {146#' example.wide147#' }148""example.wide""","R"
149"Homeostasis","MathOnco/EvoFreq","R/EvoFreq_funcs.R",".R","65726","1471","#'@import ggplot2 dplyr bezier150 151#'@title get_evofreq152#'153#' Collect information to plot frequency dynamics154#' 155#'@param size_df Dataframe in a wide format, where each row corrsponds to a single clone, and the columns are the sizes of that clone at each timepoint156#'@param clones Array containing the clone ids. The index of each clone must correspond to the same index of the row in \code{size_df} that contains the sizes of that clone over time157#'@param parents Array containing the ids of the parent of each clone in the \code{clones} array.158#'@param fill_value Array containing information that can be used to color each clone. If NULL (the default), each clone is assigned a color. If values are a clone attribute, e.g. fitness, then the colors are assigned according to those values. The user can also provide custom colors in 3 ways: 1) hexcode; 2) rgb values as a string, with each value being a the intensity of the color channel, each separated by commas, e.g. ""255, 10, 128""; 3) Any of the named in colors in R, which can be found with \code{\link[grDevices]{colors}}159#'@param fill_range Array containing the minimum and maximum values to set the range of colors. If NULL (the default), the range is determined directly from \code{fill_value}.160#'@param time_pts Array containing the name of the timepoints. If NULL, then the name of timepoints will be a sequence from 1 to the number of columns in \code{size_df}.161#'@param clone_cmap Colormap to use for the clones. For a list of available colormaps, see https://github.com/bhaskarvk/colormap.162#'@param threshold The minimum frequency of clones to be plotted. Clones with with a frequency below this value will not be plotted163#'@param scale_by_sizes_at_time Boolean defining whether or not the plot should represent the size or frequency of each clone at each timepoint. If TRUE, the sizes are scaled by the maximum size at each timepoint, and the plot thus represents the clonal frequencies at each timepoint. If FALSE, the sizes are scaled using the maximum size in \code{size_df}, thus reflecting relative population sizes164#'@param data_type String defining what kind of information is in size_df. If ""size"", then the values in \code{size_df} are the population sizes. If ""mutation"", the values are the frequencies, between 0 and 1, of each mutation in the population over time165#'@param interpolation_steps Integer defining the number of knots to use in the spline interpolation used to fill in the gaps between observed population sizes. For sparse data, this smooths out the curves in the plot. Not recommended if the data is dense, as this is slow and may not have noticable effects166#'@param interp_method String identifying the interpolation method to use. Either ""bezier"", or a method used by \code{\link[stats]{splinefun}}167#'@param fill_gaps_in_size Boolean defining whether or not missing sizes should be filled in168#'@param test_links Make sure clone does not have the same id as it's parent. If true, it can cause infinite recursion. 169#'@param add_origin Boolean defining whether or not to add origin positions to founder clones, even if not present in the data. Best for sparse observed data170#'@param tm_frac Value between 0 and 1 that determines where the maximum growth rate is in the inferred origin sizes. Lower values result in earlier maximum growth171#'@param rescale_after_thresholding Boolean determining if frequencies should be rescaled after thresholding, so that frequencies are based on what was above the threshold.172#'@param shuffle_colors Boolean determining if colors should be shuffled before being assigned to each clone. Only applies when fill_value = NULL173#'@return Formatted dataframe called a ""freq_frame"" containing the information needed to plot the frequency dynamics over time.174#'175#'@examples176#'\donttest{177#' data(""example.easy.wide"")178#' ### Split dataframe into clone info and size info using fact timepoint column names can be converted to numeric values179#' time_col_idx <- suppressWarnings(which(! is.na(as.numeric(colnames(example.easy.wide)))))180#' size_df <- example.easy.wide[, time_col_idx]181#' parents <- example.easy.wide$parents182#' clones <- example.easy.wide$clones183#' 184#' ### Default is to plot size185#' freq_frame <- get_evofreq(size_df, clones, parents)186#' evo_p_by_size <- plot_evofreq(freq_frame)187#' 188#' ### Can also plot frequency by setting scale_by_sizes_at_time = TRUE.189#' freq_frame <- get_evofreq(size_df, clones, parents, scale_by_sizes_at_time = TRUE)190#' evo_p_by_freq <- plot_evofreq(freq_frame)191#' 192#' ### Default is to mildly smooth corners, but this can be turned by setting interpolation_steps = 0193#' freq_frame <- get_evofreq(size_df, clones, parents, interpolation_steps = 0)194#' raw_evo_p <- plot_evofreq(freq_frame)195#' 196#' ### Several other methods to smooth corners, including using Bezier curves. However, Bezier curves dont represent the data as accurately as the methods that use splinefun, i.e. c(""fmm"", ""periodic"", ""natural"", ""monoH.FC"", ""hyman"")197#' freq_frame <- get_evofreq(size_df, clones, parents, interp_method = ""bezier"")198#' bez_evo_p <- plot_evofreq(freq_frame)199#' 200#' ### Data can also be provided as mutaiton frequencies by setting data_type = ""mutation""201#' mutation_count_df <- get_mutation_df(size_df, clones, parents)202#' freq_frame <- get_evofreq(mutation_count_df, clones, parents, data_type = ""mutation"")203#' evo_p_from_mutation <- plot_evofreq(freq_frame)204#' 205#' ### Input needs to be in wide format, but can be converted to long format data to wide format using \code{\link{long_to_wide_freqframe}}206#' wide_df_info <- long_to_wide_freqframe(long_pop_sizes_df = example.easy.long.sizes, time_col_name = ""Time"", clone_col_name = ""clone"", parent_col_name = ""parent"", size_col_name = ""Size"", edges_df = example.easy.long.edges)207#' clones_from_long <- wide_df_info$clones208#' parents_from_long <- wide_df_info$parents209#' size_df_from_long <- wide_df_info$wide_size_df210#' freq_frame <- get_evofreq(size_df_from_long, clones_from_long, parents_from_long)211#' evo_p_from_long <- plot_evofreq(freq_frame)212#'213#' ### Setting of colors can be done when getting the freq_frame, or by updating the color later using \code{\link{update_colors}}. For a list of available colormaps, see https://github.com/bhaskarvk/colormap.214#' ### Default colormap is rainbow_soft, but this can be changed using the \code{clone_cmap} argument. 215#' jet_freq_frame <- get_evofreq(size_df, clones, parents, clone_cmap = ""jet"")216#' jet_evo_p <- plot_evofreq(jet_freq_frame)217#' 218#' ### Can color each clone by an attribute by providing a \code{fill_value}. Default colormap is viridis, but this can be changed using the \code{clone_cmap} argument219#' fitness <- runif(length(clones))220#' fitness_freq_frame <- get_evofreq(size_df, clones, parents, fill_value = fitness)221#' fitness_evo_p <- plot_evofreq(fitness_freq_frame)222#' 223#' ### The user can also provide custom colors for each clone, which will need to be passed into the \code{fill_value} argument224#' ### Custom colors can be defined using RGB values. Each color should be a string specifying the color channel values, separated by commas.225#' rgb_clone_colors <- sapply(seq(1, length(clones)), function(x){paste(sample(0:255,size=3,replace=TRUE),collapse="","")})226#' rgb_freq_frame <- get_evofreq(size_df, clones, parents, rgb_clone_colors)227#' rgb_evo_p <- plot_evofreq(rgb_freq_frame)228#' 229#' ### Custom colors can also be any of the named colors in R. A list of the colors can be found with \code{colors()}230#' named_clone_colors <- sample(colors(), length(clones), replace = FALSE)231#' named_freq_frame <- update_colors(rgb_freq_frame, clones = clones, fill_value = named_clone_colors)232#' named_evo_p <- plot_evofreq(named_freq_frame)233#' 234#' ### Custom colors can also be specified using hexcode235#' hex_clone_colors <- c(""#614099ff"", ""#1d347eff"", ""#94558aff"", ""#c96872ff"", ""#f1884dff"", ""#e8fa5bff"", ""#042333ff"",""#f9bb41ff"")236#' hex_freq_frame <- update_colors(rgb_freq_frame, clones = clones, fill_value = hex_clone_colors)237#' hex_evo_p <- plot_evofreq(hex_freq_frame)238#'239#' ### Can revert back to original colors240#'freq_frame_default_color <- update_colors(fitness_freq_frame, clones=clones)241#'default_cmap_evo_p <- plot_evofreq(freq_frame_default_color)242#'}243#'@export244get_evofreq <- function(size_df, clones, parents, fill_value=NULL, fill_range = NULL, time_pts=NULL, clone_cmap=NULL, threshold=0.01, scale_by_sizes_at_time = FALSE, data_type=""size"", interpolation_steps = 20, interp_method = ""monoH.FC"", fill_gaps_in_size = FALSE, test_links=TRUE, add_origin=FALSE, tm_frac=0.6, rescale_after_thresholding=FALSE, shuffle_colors=FALSE){245 # # ## FOR TESTING ###246 # data(""example.easy.wide.with.attributes"")247 # ### Split dataframe into clone info and size info using fact timepoint column names can be converted to numeric values248 # time_col_idx <- suppressWarnings(which(! is.na(as.numeric(colnames(example.easy.wide)))))249 # size_df <- example.easy.wide[, time_col_idx]250 # parents <- example.easy.wide$parents251 # clones <- example.easy.wide$clones252 # fill_value <- clone_attribute_colors253 254 # parents <- clone_df$Parent255 # clones <- row.names(clone_df)256 # time_pts <- as.numeric(colnames(clone_df))257 # time_pts <- which(!is.na(time_pts))258 # size_df <- clone_df[time_pts]259 # fill_value <- clone_df$Drivers260 261 ### HAL262 # size_df <- hal_info$size_df263 # clones <- hal_info$clones264 # parents <- hal_info$parents265 # fill_value <- NULL266 267 # threshold <- 0.01268 # clone_cmap <- NULL269 # time_pts <- NULL270 # fill_range <- NULL271 # scale_by_sizes_at_time <- FALSE272 # interpolation_steps <- 10273 # fill_gaps_in_size <- FALSE274 # test_links <- TRUE275 # data_type <- ""size""276 # interp_method <- ""monoH.FC""# ""bezier""277 # add_origin <- FALSE278 # tm_frac <- 0.6279 # rescale_after_thresholding <- FALSE 280 # shuffle_colors <- FALSE281 # # # ###282 if(!is.null(fill_value)){283 fill_name <- colnames(fill_value) ### Value was passed in as a single column dataframe284 if(is.null(fill_name)){285 paresed_fill_name <- deparse(substitute(fill_value))286 print(paresed_fill_name)287 fill_name <- get_argname(paresed_fill_name)288 print(fill_name)289 }290 291 attribute_df <- data.frame(""clone_id""=clones)292 attribute_df[fill_name] <- fill_value293 }else{294 attribute_df <- NULL295 fill_name <- NULL296 }297 298 og_time_pts <- colnames(size_df)299 to_plot_df <- filter_data(size_df = size_df, clones = clones, parents = parents, time_pts = time_pts, attribute_df = attribute_df, threshold = threshold, scale_by_sizes_at_time = scale_by_sizes_at_time, data_type = data_type, fill_gaps_in_size = fill_gaps_in_size, test_links=test_links, add_origin=add_origin, tm_frac=tm_frac, rescale_after_thresholding=rescale_after_thresholding)300 clones <- to_plot_df$clones301 parents <- to_plot_df$parents302 freq_mat <- to_plot_df$freq_mat303 max_mutation_size <- to_plot_df$max_size304 time_pt_names <- to_plot_df$og_colnames305 time_pt_df <- data.frame(""x""=as.numeric(colnames(freq_mat)), ""Time_label""=time_pt_names) ### Have to make x numeric. Otherwise, float column names are converted to strings306 307 if(!is.null(attribute_df)){308 attribute_df <- to_plot_df$attributes309 }310 311 if(!is.null(fill_name)){312 fill_value <- attribute_df[,fill_name]313 if(is.null(fill_range)){314 fill_range <- range(fill_value, na.rm = TRUE)315 }316 }317 318 cat(""\n"")319 print(""Getting Plot Positions"")320 time_pts <- colnames(freq_mat)321 plot_pos_df <- get_pos(clones, parents, freq_mat, as.numeric(time_pts))322 323 324 if(interpolation_steps > 0){325 cat(""\n"")326 print(""Smoothing Polygons"")327 plot_pos_df <- smooth_pos(plot_pos_df, n_intermediate_steps = interpolation_steps, interp_method=interp_method)328 }329 plot_pos_df$extinction_time <- max(plot_pos_df$x)330 331 true_time_pt_idx <- which(colnames(freq_mat) %in% og_time_pts)332 for(cidx in seq(1, length(clones))){333 clone_freq <- freq_mat[cidx, true_time_pt_idx]334 zero_idx <- which(clone_freq==0)335 origin_time <- which(clone_freq!=0)[1]336 time_dif <- origin_time - zero_idx337 if(any(time_dif < 0)){338 extinction_time_idx <- zero_idx[which(time_dif < 0)[1]]339 # extinction_time <- as.numeric(names(time_dif)[extinction_time_idx])340 extinction_time <- og_time_pts[extinction_time_idx]341 cidx_in_pos_df <- which(plot_pos_df$clone_id == clones[cidx])342 plot_pos_df$extinction_time[cidx_in_pos_df] <- extinction_time343 }344 }345 346 ### Supply attribute name since using deparse inside get_evofreq will return fill_value for the name of the attribute347 plot_pos_df <- update_colors(evo_freq_df = plot_pos_df, clones = clones, fill_value = fill_value, clone_cmap = clone_cmap, fill_range = fill_range, fill_name=fill_name, shuffle_colors = shuffle_colors)348 349 if(!scale_by_sizes_at_time){350 plot_pos_df$y <- plot_pos_df$y*max_mutation_size351 plot_pos_df$y_label <- ""Population Size""352 353 }else{354 plot_pos_df$y_label <- ""Frequency""355 }356 357 358 plot_pos_df$row_id <- seq(1, nrow(plot_pos_df))359 updated_time_pts <- unique(plot_pos_df$x)360 ###Replace closest x value with original timepoint value before merge361 for(tpt in time_pt_df$x){362 closest_x_idx <- which.min(abs(as.numeric(tpt) - updated_time_pts))363 closest_x <- updated_time_pts[closest_x_idx]364 plot_pos_df$x[which(plot_pos_df$x==closest_x)] <- as.numeric(tpt)365 }366 367 plot_pos_df <- merge(plot_pos_df, time_pt_df, by=""x"", all=TRUE)368 plot_pos_df <- plot_pos_df[order(plot_pos_df$row_id), ]369 370 return(plot_pos_df)371}372 373scale_values <- function(x, out_range=c(0, 1)){374 a <- min(out_range)375 b <- max(out_range)376 in_min <- min(x, na.rm = TRUE)377 in_max <- max(x, na.rm = TRUE)378 379 scaled_x <- (b-a)*(x-in_min)/(in_max - in_min) + a380 381 return(scaled_x)382}383 384get_genomes <- function(c_list, p_list, out=""binary""){385 ###out=binary or bases386 n_clones <- length(c_list)387 if(out==""binary""){388 gene_mat <- matrix(0, nrow=n_clones, ncol=n_clones)389 }else{390 bases <- c(""A"", ""C"", ""T"", ""G"")391 gene_mat <- matrix(sample(bases, n_clones**2, replace = TRUE), nrow=n_clones, ncol=n_clones)392 }393 394 row.names(gene_mat) <- c_list395 for(i in seq(c_list)){396 cid <- c_list[i]397 children_idx <- get_all_idx(cid, c_list, p_list)398 if(out==""binary""){399 gene_mat[children_idx, i] <- 1400 }else{401 gene_mat[children_idx, i] <- sample(bases, 1)402 }403 }404 405 return(as.matrix(gene_mat))406}407 408traverse <- function(c, clones, parents, fnx){409 f <- match.fun(fnx)410 f(c, clones)411 for(n in get_children(c, clones, parents)){412 traverse(n, clones, parents, fnx)413 }414}415 416fill_in_gaps <- function(size_array){417 #'418 #'419 ### FOR TESTING ###420 # og_size_array <- c(0,0,0,1,2,3,4,0,0,0, 1, 2, 3, 0, 0, 1, 2, 0, 0)421 # size_array <- og_size_array422 #######423 424 425 above_zero_size_idx <- which(size_array > 0) ### indices for size array426 time_between_non_zero_idx <- diff(above_zero_size_idx) ## If there are no gaps in sizes, all these should be 1. Dif is difference between idx and idx + 1 427 left_side_of_above_0_gap_idx <- above_zero_size_idx[which(time_between_non_zero_idx > 1)] ### which indices in the original array are to the left side of a gap428 429 for(lidx in left_side_of_above_0_gap_idx){430 left_size <- size_array[lidx]431 432 temp_above_zero_idx <- which(size_array>0)433 right_non_zero_idx <- temp_above_zero_idx[which(temp_above_zero_idx>lidx)][1]434 right_size <- size_array[right_non_zero_idx]435 436 m <- (left_size-right_size)/(lidx-right_non_zero_idx)437 b <- left_size - m*lidx 438 time_btwn <- seq(lidx+1, right_non_zero_idx - 1)439 sizes_btwn <- m*time_btwn + b440 size_array[time_btwn] <- sizes_btwn441 }442 return(size_array)443}444 445check_for_missing_links_fnx <- function(clones, parents){446 for(cid in clones){447 get_ancestor_idx_evo_freq(cid, clone_list = clones, parent_list = parents)448 }449 return(""No missing links"")450}451 452check_and_update_edges <- function(clones, parents, check_for_missing_links=TRUE){453 ##### FOR TESTING ###454 # clones <- clone_list455 # parents <- parent_list456 ######457 458 ### Make sure each clone only appears once in edge list459 if(all(as.numeric(table(clones)) != 1)){460 warning(""Clone occurs more than once in edgelist"")461 }462 463 ### Update so that parent of root has an id not in the clone_list. 464 clone_as_parent_idx <- which(clones==parents)465 if(length(clone_as_parent_idx)>0){466 if(length(clone_as_parent_idx)>1){467 warning(""More than 1 clone is its own parent. Cannot determine root"")468 }469 rand_parent_id <- stats::runif(1, min = -1, max=0)470 parents[clone_as_parent_idx] <- rand_parent_id471 }472 473 ### Determine if root can be found ###474 root_id <- get_root_id(clones, parents)475 476 return(list(""updated_clones""=clones, ""updated_parents""=parents))477}478 479check_freq_mat <- function(freq_mat, clones, parents){480 ### FOR TESTING ###481 # freq_mat <- mutation_df482 # clones <- clone_list483 # parents <- parent_list484 #########485 486 ### Make sure total number children carrying the mutation ###487 for(i in seq(1, length(clones))){488 cln <- clones[i]489 cln_freq <- freq_mat[i, ]490 all_children_idx <- get_all_idx(cln, clones, parents)491 children <- clones[all_children_idx]492 children_idx <- all_children_idx[clones[all_children_idx] != cln]493 if(length(children == 0)){494 next495 }496 children_df <- freq_mat[children_idx, ]497 parent_greater_than_children <- apply(children_df, 1, function(x){all(cln_freq>=x)})498 if(all(parent_greater_than_children)==FALSE){499 too_big_idx <- which(parent_greater_than_children==FALSE)500 clones_too_big <- clones[too_big_idx]501 stop(paste(cln, ""has descendent mutations greater than that are greater than it's size. Descendents are:"", clones_too_big))502 }503 }504}505 506#'\code{get_mutation_df} Converts sizes to the frequency of the mutations in the population507#'@inheritParams get_evofreq508#'@param clone_size_df Dataframe in a wide format, where each row corrsponds to a single clone, and the columns are the sizes of that clone at each timepoint509#'@examples510#'data(""example.easy.wide"")511#'### Split dataframe into clone info and size info using fact timepoint column names can be converted to numeric values512#'time_col_idx <- suppressWarnings(which(! is.na(as.numeric(colnames(example.easy.wide)))))513#'attribute_col_idx <-suppressWarnings(which(is.na(as.numeric(colnames(example.easy.wide)))))514#'attribute_df <- example.easy.wide[, attribute_col_idx]515#'size_df <- example.easy.wide[, time_col_idx]516#'parents <- example.easy.wide$parents517#'clones <- example.easy.wide$clones518#'### Size data can be converted to mutation counts for additional analyses, like calculating mutation frequencies519#'mutation_count_df <- get_mutation_df(size_df, clones, parents)520#'total_chromosomes <- 2*colSums(size_df) ### Calculation assumes that mutations are on a single chromosome521#'allele_frequency_df <- sweep(mutation_count_df,2,total_chromosomes,""/"")522#'@export523get_mutation_df <- function(clone_size_df, clones, parents){524 ### FOR TESTING ####525 # clone_size_df <- as.data.frame(size_df)526 # clones <- clone_list527 # parents <- parent_list528 ####529 cat(""\n"")530 print(""Getting Mutation Counts"")531 greater_than_one_time_pt <- ncol(clone_size_df)>1532 mutation_df <- matrix(NA, nrow = nrow(clone_size_df), ncol = ncol(clone_size_df))533 534 pb <- utils::txtProgressBar(min = 1, max = nrow(clone_size_df), style = 3)535 for(i in seq(1, length(clones))){536 utils::setTxtProgressBar(pb, i)537 cln <- clones[i]538 all_children_idx <- get_all_idx(cln, clones, parents)539 540 child_sizes_df <- clone_size_df[all_children_idx, ]541 if(greater_than_one_time_pt){542 total_cells_carrying_mutation <- colSums(child_sizes_df) ### recursion already includes adding the current clone 543 }else{544 total_cells_carrying_mutation <- sum(child_sizes_df)545 }546 547 mutation_df[i, ] <- total_cells_carrying_mutation548 549 }550 551 # check_freq_mat(mutation_df, clones, parents)552 colnames(mutation_df) <- colnames(clone_size_df)553 rownames(mutation_df) <- rownames(clone_size_df)554 return(mutation_df)555}556 557predict_sizes <- function(x1, y1, x2, y2, new_x){558 ### FOR Testing ###559 # y1 <- 1.00560 # y2 <- 0.02561 # x1 <- 0562 # x2 <- 30563 # new_x <- 15564 ####565 566 567 fit <- stats::lm(c(y1, y2) ~ c(x1, x2))568 new_y <- fit$coefficients[[1]] + fit$coefficients[[2]]*new_x569 return(new_y)570}571 572make_pointy <- function(dull_mat){573 ### FOR TESTING ###574 # dull_mat <- new_mat575 ####576 ### Make that polygons do not have abrupt starts and stops577 ### Set size to 1 for intermediate timestep before origin578 ### If clone went extinct, set intemediate time step after extinction to 1579 580 581 time_pts <- colnames(dull_mat)582 max_time <- max(time_pts)583 min_time <- min(time_pts)584 for(cidx in seq(nrow(dull_mat))){585 # cidx <- 2586 clone_origin_time_idx <- which(dull_mat[cidx, ] > 0)[1]587 clone_origin_time <- time_pts[clone_origin_time_idx]588 if(clone_origin_time > min_time){589 dull_mat[cidx, clone_origin_time_idx-1] <- 0.0001590 }591 time_zero_idx <- which(dull_mat[cidx, ] == 0)592 if(length(time_zero_idx) > 0){593 594 if(any(time_zero_idx > clone_origin_time_idx)){595 extinction_time <- names(which(time_zero_idx > clone_origin_time_idx)[1])596 extinction_time_idx <- which(time_pts == extinction_time)597 dull_mat[cidx, extinction_time_idx] <- 0.0001 598 599 } 600 }601 }602 return(dull_mat)603}604 605interp_mut_mat <- function(mut_mat, time_pts= NULL){606 #For each clone, estimate intermediate sizes607 ### FOR TESTING ###608 # mut_mat <- og_mut_mat609 # time_pts <- og_time_pts610 ####611 612 if(is.null(time_pts)){613 time_pts <- seq(1, ncol(mut_mat))614 }615 616 intermediate_vals <- mapply(i=seq(2, length(time_pts)), FUN=function(i){(time_pts[i]-time_pts[i-1])/2 + time_pts[i-1]})617 new_times <- rep(NA, length(time_pts)+length(intermediate_vals))618 og_time_idx <- seq(1, length(new_times),2)619 new_time_idx <- seq(2, length(new_times)-1,2)620 621 new_times[og_time_idx] <- time_pts622 new_times[new_time_idx] <- intermediate_vals623 n_new_times <- length(new_times)624 ## Increase width of matrix625 626 new_mat <- matrix(NA, nrow=nrow(mut_mat), ncol=n_new_times)627 colnames(new_mat) <- new_times628 row.names(new_mat) <- row.names(mut_mat)629 630 # max_time <- max(time_pts)631 # min_time <- min(time_pts)632 633 for(cidx in seq(nrow(new_mat))){634 # cidx <- 2635 clone_sizes <-as.numeric( mut_mat[cidx, ])636 intermediate_sizes <- mapply(i=seq(2, length(time_pts)), function(i){637 # i <- 2638 if(clone_sizes[i]== 0 | clone_sizes[i-1]== 0){639 return(0)640 }641 642 #643 fit <- stats::lm(c(clone_sizes[i], clone_sizes[i-1]) ~ c(time_pts[i], time_pts[i-1]))644 new_y <- fit$coefficients[[1]] + fit$coefficients[[2]]* intermediate_vals[i-1]645 # new_y <- fit$coefficients[[1]] + fit$coefficients[[2]]* intermediate_vals[i]646 return(new_y)647 })648 # intermediate_sizes <- clone_sizes[1:length(clone_sizes)-1] + diff(clone_sizes)/2649 new_mat[cidx, og_time_idx] <- clone_sizes650 new_mat[cidx, new_time_idx] <- intermediate_sizes651 652 }653 654 new_mat <- make_pointy(new_mat)655 return(new_mat)656 657}658 659order_clones_at_time <- function(clones_at_time, parents, clones){660 ### possible that parent and child arose within the same period. Need to ensure that parents are drawn first661 662 #### For testing ###663 # clones_at_time <- parents_at_time664 # parents <- parents665 # clones <- clones666 ####667 root_id <- get_root_id(clones, parents)668 n_clones <- length(clones_at_time)669 n_ancestors_in_t <- rep(NA, n_clones)670 671 for(ct_i in seq(n_clones)){672 cid <- clones_at_time[ct_i]673 if(cid==root_id){674 ### makes sure root is always drawn first675 n_ancestors_in_t[ct_i] <- -1676 next()677 }678 ancestor_idx <- get_ancestor_idx_evo_freq(cid, clones, parents)679 if(is.null(ancestor_idx)){680 n_ancestors_in_t[ct_i] <- 0681 }else{682 n_ancestors_in_t[ct_i] <- length(ancestor_idx)683 }684 }685 new_order <- order(n_ancestors_in_t, decreasing = FALSE)686 clones_at_time_ordered <- clones_at_time[new_order]687 688 return(clones_at_time_ordered)689}690 691distribute_x <- function(x){692 #### For testing ####693 # x <- seq(1, 6)694 # x <- c(1)695 # x <- c(1800, 2008, 2166, 2202, 3564, 3694)696 #####################697 698 ### Put even numbers on one side of first point, odd numbers on the other699 x_idx <- seq(2, length(x))700 even_idx <- x_idx[x_idx %% 2 == 0]701 odd_idx <- x_idx[x_idx %% 2 != 0]702 new_order_idx <- c(rev(even_idx), 1, odd_idx)703 new_order <- x[new_order_idx]704 return(new_order)705}706 707filter_freq_mat <- function(clones, parents, freq_mat, threshold, attr_vals=NULL){708 ### FOR TESTING ###709 # attr_vals <- attribute_vals710 ##711 max_freqs <- apply(freq_mat, 1, max)712 above_thresh_idx <- which(max_freqs >= threshold)713 714 filtered_freq_mat <- freq_mat[above_thresh_idx, ]715 filtered_clones <- clones[above_thresh_idx]716 filtered_parents <- parents[above_thresh_idx]717 if(!is.null(attr_vals)){718 filtered_attr_vals <- attr_vals[above_thresh_idx]719 }else{720 filtered_attr_vals <- NULL721 }722 723 return_vals <- list('freq_mat'=filtered_freq_mat, 'clones'= filtered_clones, 'parents'=filtered_parents, 'attr_vals'=filtered_attr_vals, 'thresh_idx'=above_thresh_idx)724}725 726linear_interp_point <- function(x1, x2, y1, y2){727 #### FOR TESTING ####728 # x1 <- new_x[fp_idx-1]729 # x2 <- new_x[fp_idx+1]730 # y1 <- new_forward[fp_idx-1]731 # y2 <- new_forward[fp_idx+1]732 #####733 m <- (y1-y2)/(x1-x2)734 b <- y1 - m*x1735 middle_x <- 0.5*(x1 + x2)736 new_y <- m*middle_x + b737 return(new_y)738}739 740get_linear_m_and_b <- function(x1, x2, y1, y2){741 #### FOR TESTING ####742 # x1 <- new_x[fp_idx-1]743 # x2 <- new_x[fp_idx+1]744 # y1 <- new_forward[fp_idx-1]745 # y2 <- new_forward[fp_idx+1]746 #####747 m <- (y1-y2)/(x1-x2)748 b <- y1 - m*x1749 return(list(""m""=m, ""b""=b))750}751 752smooth_pos <- function(sparse_pos_df, n_intermediate_steps=20, interp_method = ""monoH.FC""){753 ###FOR TESTING ##754 # sparse_pos_df <- evo_pos755 # n_intermediate_steps <- 3756 # interp_method <- ""bezier""757 ###758 759 interp_df_list <- list()760 unique_clones <- unique(sparse_pos_df$clone_id)761 n_clones <- length(unique_clones)762 763 all_x <- unique(sparse_pos_df$x)764 x_range <- range(all_x, na.rm = TRUE)765 all_new_x <- seq(x_range[1], x_range[2], length.out = length(all_x)*n_intermediate_steps)766 767 pb <- utils::txtProgressBar(min = 0, max = n_clones, style = 3)768 for(i in seq(n_clones)){769 # print(i)770 # i <- 4771 cid <- unique_clones[i]772 utils::setTxtProgressBar(pb, i)773 clone_pos_df <- subset(sparse_pos_df, clone_id==cid)774 nx <- length(unique(clone_pos_df$x))775 776 if(nx >= 4){777 778 forward_idx <- which(duplicated(clone_pos_df$x)==FALSE)779 forward_df <- clone_pos_df[forward_idx,]780 ###GO BACK TO THIS??781 # new_x <- seq(min(forward_df$x), max(forward_df$x), length.out = length(forward_df$x)*n_intermediate_steps)782 ###783 784 clone_x_range <- range(forward_df$x)785 new_x <- all_new_x[all_new_x >= clone_x_range[1] & all_new_x <= clone_x_range[2]]786 787 ### need at least 4 points for spline interpolation788 # print(interp_method)789 # print(interp_method %in% c(""fmm"", ""periodic"", ""natural"", ""monoH.FC"", ""hyman""))790 if(interp_method %in% c(""fmm"", ""periodic"", ""natural"", ""monoH.FC"", ""hyman"")){791 # print(""FW"")792 # interp_method <- ""monoH.FC""793 func = stats::splinefun(x=forward_df$x, y=forward_df$y, method=interp_method, ties = mean)794 new_forward <- func(new_x)795 new_forward_x <- new_x796 797 }else if(interp_method==""bezier""){798 bezier_forward_x <- seq(0, 1, length.out = length(new_x))799 forward_bezier_points <- bezier::bezier(t=bezier_forward_x, p=forward_df[,c(""x"",""y"")])800 new_forward <- forward_bezier_points[, 2]801 new_forward_x <- forward_bezier_points[, 1]802 803 }else if(interp_method==""bezier_curve_fit""){804 bezier_forward_x <- seq(0, 1, length.out = length(new_x))805 bz_forward_fit <- suppressWarnings(bezier::bezierCurveFit(as.matrix(forward_df[,c(""x"",""y"")]), na.fill = TRUE, maxiter=500, max.rse.percent.change = 0.2, fix.start.end = TRUE))806 forward_bezier_points <- bezier::bezier(t=bezier_forward_x, p=bz_forward_fit$p)807 new_forward <- forward_bezier_points[, 2]808 new_forward_x <- forward_bezier_points[, 1]809 }else if(interp_method==""loess""){810 fit <- loess(forward_df$y ~ forward_df$x)811 new_forward <- predict(fit, new_x)812 new_forward_x <- new_x813 }814 815 816 reverse_idx <- which(duplicated(clone_pos_df$x)==TRUE)817 reverse_df <- clone_pos_df[reverse_idx,]818 819 if(interp_method %in% c(""fmm"", ""periodic"", ""natural"", ""monoH.FC"", ""hyman"")){820 # print(""REVERSE"")821 rev_func = stats::splinefun(x=rev(reverse_df$x), y=rev(reverse_df$y), method=interp_method, ties = mean) ### x has to increase822 new_reverse <- rev_func(new_x)823 new_reverse_x <- rev(new_x)824 825 }else if(interp_method==""bezier""){826 bezier_reverse_x <- seq(0, 1, length.out = length(new_x))827 reverse_bezier_points <- bezier::bezier(t=bezier_reverse_x, p= reverse_df[seq(nrow(reverse_df), 1),c(""x"",""y"")])828 new_reverse <- reverse_bezier_points[, 2]829 new_reverse_x <- rev(reverse_bezier_points[, 1])830 }else if(interp_method==""bezier_curve_fit""){831 bezier_reverse_x <- seq(0, 1, length.out = length(new_x))832 bz_reverse_fit <- suppressWarnings(bezier::bezierCurveFit(as.matrix(reverse_df[seq(nrow(reverse_df), 1),c(""x"",""y"")]), na.fill = TRUE, maxiter=500, max.rse.percent.change = 0.2, fix.start.end = TRUE))833 reverse_bezier_points <- bezier::bezier(t=bezier_reverse_x, p=bz_reverse_fit$p)834 new_reverse <- reverse_bezier_points[, 2]835 new_reverse_x <- rev(reverse_bezier_points[, 1]) 836 }else if(interp_method==""loess""){837 rev_fit <- loess(rev(reverse_df$y) ~ rev(reverse_df$x))838 new_reverse <- predict(rev_fit, new_x)839 new_forward_x <- rev(new_x)840 }841 842 843 #### Make sure that there are no times where the bottom and top cross. Bottom (forward) should always be less than top (reverse)844 rev_forward_dist <- new_reverse - new_forward845 flip_idx <- which(rev_forward_dist < 0)846 n_new_x <- length(new_x)847 if(length(flip_idx) > 0){848 # fp_idx <- 2849 for(fp_idx in flip_idx){850 updated_dist <- new_reverse[fp_idx] - new_forward[fp_idx]851 if(updated_dist > 0){852 #### position was already updated853 next()854 }855 if(fp_idx > n_new_x | n_new_x <= 1 | fp_idx == 1){856 ### first and last point are the new pointy ends. This just flips them to make a triangle857 replace_forward <- new_reverse[fp_idx]858 replace_reverse <- new_forward[fp_idx]859 new_forward[fp_idx] <- replace_forward860 new_reverse[fp_idx] <- replace_reverse861 862 }else{863 864 865 ### For forward: get next positive difference (should have been negative because forward should be < reverse). Connect to that point. Will also need to update any other points along that path866 next_seq <- seq(fp_idx+1, n_new_x)867 # next_positive_idx <- fp_idx + which(rev_forward_dist[fp_idx+1:n_new_x] > 0)[1]868 next_positive_idx <- fp_idx + which(rev_forward_dist[next_seq] > 0)[1]869 if(is.na(next_positive_idx)){870 ###ALL positions after this one need to be flipped871 replace_forward <- new_reverse[next_seq]872 replace_reverse <- new_forward[next_seq]873 new_forward[next_seq] <- replace_forward874 new_reverse[next_seq] <- replace_reverse875 break()876 }877 replace_seq <- seq(fp_idx, next_positive_idx - 1)878 879 forward_lm <- get_linear_m_and_b(new_x[fp_idx-1], new_x[next_positive_idx], new_forward[fp_idx-1], new_forward[next_positive_idx])880 new_forward[replace_seq] <- forward_lm$m*new_x[replace_seq] + forward_lm$b881 # new_forward[replace_seq] <- func(forward_lm$m*new_x[replace_seq] + forward_lm$b)882 883 ### For reverse: get next negative difference. Connect to that point. Will also need to update any other points along that path884 reverse_lm <- get_linear_m_and_b(new_x[fp_idx-1], new_x[next_positive_idx], new_reverse[fp_idx-1], new_reverse[next_positive_idx])885 new_reverse[replace_seq] <- reverse_lm$m*new_x[replace_seq] + reverse_lm$b886 887 # replace_forward <- linear_interp_point(new_x[fp_idx-1], new_x[fp_idx+1], new_forward[fp_idx-1], new_forward[fp_idx+1])888 # replace_reverse <- linear_interp_point(new_x[fp_idx-1], new_x[fp_idx+1], new_reverse[fp_idx-1], new_reverse[fp_idx+1])889 }890 891 # if(replace_reverse - replace_forward ){892 # new_forward[fp_idx] <- replace_reverse893 # new_reverse[fp_idx] <- replace_forward894 # }else{895 896 }897 }898 899 # print(c(length(new_x), length(new_forward), length(new_reverse)))900 new_forward_df <- data.frame('x'=new_x, 'y'=new_forward, 'clone_id'=clone_pos_df$clone_id[1], 'parent'=clone_pos_df$parent[1], 'origin_time'=clone_pos_df$origin_time[1], 'draw_order'=clone_pos_df$draw_order[1])901 new_reverse_df <- data.frame('x'=rev(new_x), 'y'=rev(new_reverse), 'clone_id'=clone_pos_df$clone_id[1], 'parent'=clone_pos_df$parent[1], 'origin_time'=clone_pos_df$origin_time[1], 'draw_order'=clone_pos_df$draw_order[1])902 903 # new_forward_df <- data.frame('x'=new_forward_x, 'y'=new_forward, 'clone_id'=clone_pos_df$clone_id[1], 'parent'=clone_pos_df$parent[1], 'origin_time'=clone_pos_df$origin_time[1], 'draw_order'=clone_pos_df$draw_order[1])904 # new_reverse_df <- data.frame('x'=new_reverse_x, 'y'=rev(new_reverse), 'clone_id'=clone_pos_df$clone_id[1], 'parent'=clone_pos_df$parent[1], 'origin_time'=clone_pos_df$origin_time[1], 'draw_order'=clone_pos_df$draw_order[1])905 906 new_clone_pos_df <- rbind(new_forward_df, new_reverse_df)907 908 # gp <- ggplot2::ggplot(new_clone_pos_df, ggplot2::aes(x=x, y=y)) +909 # ggplot2::geom_polygon()910 # print(gp)911 912 913 }else{914 new_clone_pos_df <- clone_pos_df[colnames(clone_pos_df) != ""size""]915 }916 interp_df_list[[as.character(cid)]] <- new_clone_pos_df917 }918 interp_df <- do.call(rbind, interp_df_list)919 920 921 ### interpolation may cause points to go out of bounds922 interp_df$y[interp_df$y < 0] <- 0923 interp_df$y[interp_df$y > 1] <- 1924 925 926 return(interp_df)927 928}929 930#'\code{plot_evofreq} Plots the frequency dynamics using ggplot. Can also use info to create animations 931#'@param freq_frame Proper dataframe returned by \code{\link{get_evofreq}}932#'@param n_time_pts Integer defining how many time points to plot, evenly spaced between \code{start_time} and \code{end_time}. If NULL, all timepoints are plotted933#'@param start_time Integer defining the timepoint at which to start plotting frequency dynamics. If NULL, plotting will begin at the first timepoint934#'@param end_time Integer defining the timepoint at which to stop plotting frequency dynamics. If NULL, plotting will end at the last timepoint935#'@param bw width of lines surrounding polygons936#'@param bc color of lines surrounding polygons937#'@param show_axes Whether to display axes938#'@param fill_range Array containing the minimum and maximum values to set the range of colors. If NULL (the default), the range is determined directly from \code{fill_value}.939#'@return ggplot of the frequency dynamics940#'@examples941#' 942#' data(""example.easy.wide.with.attributes"")943#' ### Split dataframe into clone info and size info using fact timepoint column names can be converted to numeric values944#' time_col_idx <- suppressWarnings(which(! is.na(as.numeric(colnames(example.easy.wide.with.attributes)))))945#' attribute_col_idx <- suppressWarnings(which(is.na(as.numeric(colnames(example.easy.wide.with.attributes)))))946#' size_df <- example.easy.wide.with.attributes[, time_col_idx]947#' parents <- example.easy.wide.with.attributes$parent948#' clones <- example.easy.wide.with.attributes$clone949#' fitness <- example.easy.wide.with.attributes$fitness950#' 951#' #' 952#' ### Setting of colors can be done when getting the freq_frame, or by updating the color later using \code{\link{update_colors}}. For a list of available colormaps, see https://github.com/bhaskarvk/colormap.953#' ### Default colormap is rainbow_soft, but this can be changed using the \code{clone_cmap} argument. 954#' freq_frame <- get_evofreq(size_df, clones, parents)955#' evo_p <- plot_evofreq(freq_frame)956#' 957#' ### Can color each clone by an attribute by providing a \code{fill_value}. Default colormap is viridis, but this can be changed using the \code{clone_cmap} argument958#' fitness <- runif(length(clones))959#' fitness_freq_frame <- get_evofreq(size_df, clones, parents, fill_value = fitness)960#' fitness_evo_p <- plot_evofreq(fitness_freq_frame)961#' #' 962#' ### The user can also provide custom colors for each clone, which will need to be passed into the \code{fill_value} argument963#' ### Custom colors can be defined using RGB values. Each color should be a string specifying the color channel values, separated by commas.964#' rgb_clone_colors <- sapply(seq(1, length(clones)), function(x){paste(sample(0:255,size=3,replace=TRUE),collapse="","")})965#' rgb_freq_frame <- get_evofreq(size_df, clones, parents, rgb_clone_colors)966#' rgb_evo_p <- plot_evofreq(rgb_freq_frame)967#' 968#' ### Custom colors can also be any of the named colors in R. A list of the colors can be found with \code{colors()}969#' named_clone_colors <- sample(colors(), length(clones), replace = FALSE)970#' named_freq_frame <- update_colors(rgb_freq_frame, clones = clones, fill_value = named_clone_colors)971#' named_evo_p <- plot_evofreq(named_freq_frame)972#' 973#' ### Custom colors can also be specified using hexcode974#' hex_clone_colors <- c(""#614099ff"", ""#1d347eff"", ""#94558aff"", ""#c96872ff"", ""#f1884dff"", ""#e8fa5bff"", ""#042333ff"",""#f9bb41ff"")975#' hex_freq_frame <- update_colors(rgb_freq_frame, clones = clones, fill_value = hex_clone_colors)976#' hex_evo_p <- plot_evofreq(hex_freq_frame)977#'978#' ### Can revert back to original colors979#'freq_frame_default_color <- update_colors(fitness_freq_frame, clones=clones)980#'default_cmap_evo_p <- plot_evofreq(freq_frame_default_color)981#'982#'### Can add gganimate objects to evo_p_by_size to create the animation983#'984#'\donttest{985#' library(gganimate)986#' movie_p <- fitness_evo_p +987#' transition_reveal(x, range=range(fitness_freq_frame$x)) +988#' view_follow()989#' # print(movie_p)990#' anim_save(""evofreq_movie.gif"", movie_p)991#' }992#'@export993plot_evofreq <- function(freq_frame, n_time_pts=NULL, start_time=NULL, end_time=NULL, bw=0.05, bc=""grey75"", show_axes=TRUE, fill_range=NULL){994 995 ### FOR TESTING ###996 # n_time_pts <- NULL997 # start_time <- NULL998 # end_time <- NULL999 # bw <- 0.051000 # bc <- ""grey75""1001 # freq_frame <- hal_plot_df1002 #####1003 1004 unique_time_pts <- unique(freq_frame$x)1005 if(is.null(n_time_pts)){1006 n_time_pts <- length(unique_time_pts)1007 }1008 if(is.null(start_time)){1009 start_time <- 01010 }1011 if(is.null(end_time)){1012 end_time <- max(freq_frame$x)1013 }1014 1015 1016 # end_time <- 501017 view_df <- subset(freq_frame, x <= end_time)1018 1019 color_attribute_name <- unique(view_df$efp_color_attribute)1020 if(is.na(color_attribute_name)){1021 color_attribute_name <- ""plot_color""1022 }else{1023 if(is.null(fill_range)){1024 fill_range <- range(view_df[, color_attribute_name], na.rm = TRUE) 1025 }1026 1027 # color_df <- view_df[duplicated(view_df$plot_color)==FALSE, ]1028 # color_df <- color_df[order(color_df[color_attribute_name]), ]1029 # colorbar_colors <- color_df$plot_color1030 # 1031 }1032 1033 1034 # all_x_maxs_at_time <- aggregate(x ~ clone_id, data=subset(view_df, extinction_time >= end_time), FUN=max)1035 # max_x <- min(all_x_maxs_at_time$x) 1036 # view_df$x <- as.numeric(view_df$x)1037 time_pts_df <- unique(view_df[complete.cases(view_df[""Time_label""]), ][c(""x"", ""Time_label"")])1038 time_pts_df$x <- as.numeric(as.character(time_pts_df$x))1039 time_pts_df <- time_pts_df[order(time_pts_df$x), ]1040 time_pts_df <- subset(time_pts_df, x > 0)1041 time_pts_df$Time_label <- as.character(time_pts_df$Time_label)1042 # time_pts_df$Time_label <- as.character(as.numeric(time_pts_df$Time_label))1043 1044 # time_pts_df$Time_label[1] <- ""0"" ### TODO Why is this being set to 0?1045 1046 y_label <- unique(freq_frame$y_label)1047 1048 ggevodyn <- ggplot2::ggplot(view_df, ggplot2::aes_string(x=""x"", y=""y"", group=""draw_order"", fill=color_attribute_name)) +1049 ggplot2::geom_polygon(size=bw, color=bc) +1050 ggplot2::ylab(y_label) +1051 ggplot2::xlab(""Time"") +1052 ggplot2::theme_classic() +1053 ggplot2::scale_x_continuous(breaks=as.numeric(as.character(time_pts_df$x)), labels=as.character(time_pts_df$Time_label)) +1054 ggplot2::theme(axis.text.x = ggplot2::element_text(angle = 45, hjust = 1))1055 1056 if(color_attribute_name==""plot_color""){1057 ggevodyn <- ggevodyn + ggplot2::scale_fill_identity()1058 }else{1059 colormap_name <- unique(view_df$cmap)1060 colorbar_colors <- get_colors(100, colormap_name)1061 ggevodyn <- ggevodyn + ggplot2::scale_fill_gradientn(colours = colorbar_colors, limits=fill_range) #colormap::scale_fill_colormap(color_attribute_name, colormap=colormap_name)1062 # ggevodyn <- ggevodyn + colormap::scale_fill_colormap(color_attribute_name, colormap=colormap_name, limits=fill_range)1063 }1064 1065 if(! show_axes){1066 ggevodyn <- ggevodyn + ggplot2::theme_void()1067 }1068 return(ggevodyn)1069}1070 1071#'Convert data in long format to wide format1072#'\code{long_to_wide_freq_ready} Converts long data frame to wide format used by evofreq.1073#'@param edges_df Dataframe with 2 columns, where each row defines relation between each clone and their parent. One column contains the clone ids, and the other contains the id of that clone's parent.1074#'@param long_pop_sizes_df Dataframe with 3 columns, one containing the clone id, one containing the time points, and one containing the size of the clone at that time point. There must be a size for each clone at every time point, so it did not exist, it's size is 0.1075#'@param time_col_name String that specifies which column in long_pop_sizes_df contains the timepoint information1076#'@param clone_col_name String that specifies which column in in edges_df contains the names of the clones. I.e. the ""to"" node in a graph1077#'@param parent_col_name String that specifies which column in in edges_df contains the names of the parents of the descendent clone in the edges_df with colname clone_col_name. #I.e. the ""from"" node in a graph1078#'@param size_col_name String that specifies which column in long_pop_sizes_df contains the size of each clone for the corresponding timepoint and clone in long_pop_sizes_df1079#'@param fill_gaps_in_size Boolean specificing if gaps in sizes over time should be filled in, assuming linear changes. 1080#'@return List containing: wide_size_df, a dataframe of the sizes in wide format, where each row contains the sizes over time for a single clone; parents, a vector containing the parent id for each clone in clones; clones, a vector of clones ids correspond to each row in wide_size_df1081#'@examples1082# ### Input needs to be in wide format, but can convert long format data to wide format using \code{long_to_wide_freq_ready}1083# wide_df_info <- long_to_wide_freq_ready(long_pop_sizes_df = example.easy.long.sizes, time_col_name = ""Time"", clone_col_name = ""clone"", parent_col_name = ""parent"", size_col_name = ""Size"", edges_df = example.easy.long.edges)1084# clones_from_long <- wide_df_info$clones1085# parents_from_long <- wide_df_info$parents1086# size_df_from_long <- wide_df_info$wide_size_df1087# pos_from_long_df <- get_evofreq(size_df_from_long, clones_from_long, parents_from_long)1088# evo_p_from_long <- plot_evofreq(pos_from_long_df)1089#'@export1090long_to_wide_freq_ready <- function(edges_df, long_pop_sizes_df, time_col_name, clone_col_name, parent_col_name, size_col_name, fill_gaps_in_size=FALSE){1091 long_pop_sizes_df <- as.data.frame(long_pop_sizes_df)1092 unique_times <- unique(long_pop_sizes_df[,time_col_name])1093 unique_times <- unique_times[order(unique_times)]1094 n_time_pts <- length(unique_times)1095 1096 clones_in_edge_df <- edges_df[, clone_col_name]1097 parents_in_edge_df <- edges_df[, parent_col_name]1098 clones_in_size_df <- unique(long_pop_sizes_df[,clone_col_name])1099 clones_not_in_clone_list <- clones_in_size_df[!clones_in_size_df %in% clones_in_edge_df]1100 1101 if(length(clones_not_in_clone_list)>0){1102 if(length(clones_not_in_clone_list)>1){1103 stop(""More thant 1 clone in size data frame but not in clone list"")1104 }1105 warning(""Clone in size data frame but not in clone list. Assuming this is the root"")1106 rand_parent_id <- stats::runif(1)1107 clones_in_edge_df <- c(clones_in_edge_df, clones_not_in_clone_list)1108 parents_in_edge_df <- c(parents_in_edge_df, rand_parent_id)1109 }else{1110 # root_id <- get_root_id(clones_in_edge_df, parents_in_edge_df)1111 clone_as_parent_idx <- which(parents_in_edge_df==clones_in_edge_df)1112 if(length(clone_as_parent_idx)>0){1113 rand_parent_id <- stats::runif(1)1114 parents_in_edge_df[clone_as_parent_idx] <- rand_parent_id1115 # clones_in_edge_df <- c(clones_in_edge_df, clones_not_in_clone_list)1116 # parents_in_edge_df <- c(parents_in_edge_df, rand_parent_id)1117 if(length(clone_as_parent_idx)>1){1118 warning(""More than 1 clone is its own parent. Cannot determine root"")1119 } 1120 }1121 }1122 1123 edges_df <- data.frame(parents_in_edge_df, clones_in_edge_df)1124 names(edges_df) <- c(parent_col_name, clone_col_name)1125 1126 n_clones <- length(clones_in_edge_df)1127 size_mat <- matrix(0, nrow=n_clones, ncol=length(unique_times))1128 colnames(size_mat) <- unique_times1129 1130 size_mat_rownames <- rep(NA, n_clones)1131 cat(""\n"")1132 print(""Converting From Long to Wide Format"")1133 pb <- utils::txtProgressBar(min = 1, max = n_clones, style = 3)1134 for(i in seq(n_clones)){1135 utils::setTxtProgressBar(pb, i)1136 cid <- clones_in_edge_df[i]1137 cidx <- which(long_pop_sizes_df[clone_col_name]==cid)1138 if(length(cidx)== 0){1139 size_mat[i, ] <- 01140 }else{1141 cid_sizes <- long_pop_sizes_df[cidx, size_col_name] 1142 cid_times <- long_pop_sizes_df[cidx, time_col_name]1143 c_time_idx <- sapply(cid_times, function(x){which(unique_times==x)})1144 1145 size_mat[i, c_time_idx] <- cid_sizes1146 }1147 size_mat_rownames[i] <- cid1148 }1149 rownames(size_mat) <- size_mat_rownames1150 1151 if(fill_gaps_in_size){1152 wide_size_df <- as.data.frame(t(apply(size_mat, 1, fill_in_gaps)))1153 }else{1154 wide_size_df <- as.data.frame(size_mat)1155 }1156 1157 return(list(""wide_size_df""=wide_size_df, ""parents""=edges_df[, parent_col_name], ""clones""=edges_df[, clone_col_name]))1158 1159 1160 1161 ###Reorder edges so that they are the same as 1162 # long_pop_sizes_df <- long_pop_sizes_df[order(long_pop_sizes_df[,time_col_name]),]1163 # wide_size_df <- reshape(as.data.frame(long_pop_sizes_df), timevar = time_col_name, idvar = clone_col_name, direction = ""wide"")1164 1165 # parents <- edges_df[, parent_col_name]1166 # clones <- edges_df[, clone_col_name]1167 # clones_in_size_df <- unique(long_pop_sizes_df[,clone_col_name])1168 # clones_not_in_clone_list <- clones_in_size_df[!clones_in_size_df %in% clones]1169 # if(length(clones_not_in_clone_list)>0){1170 # if(length(clones_not_in_clone_list)>1){1171 # stop(""More thant 1 clone in size data frame but not in clone list"")1172 # }1173 # warning(""Clone in size data frame but not in clone list. Assuming this is the root"")1174 # rand_parent_id <- runif(1)1175 # clones <- c(clones, clones_not_in_clone_list)1176 # parents <- c(parents, rand_parent_id)1177 # edges_df <- data.frame(clones, parents)1178 # names(edges_df) <- c(clone_col_name, parent_col_name)1179 # }else{1180 # 1181 # }1182 1183 # wide_size_df <- merge(wide_size_df, edges_df)1184 1185 # ordered_clones <- wide_size_df[,clone_col_name]1186 # ordered_parents <- wide_size_df[,parent_col_name]1187 # ordered_size_df <- wide_size_df[,!colnames(wide_size_df) %in% c(clone_col_name, parent_col_name)]1188 # return(list(""wide_size_df""=ordered_size_df, ""parents""=ordered_parents, ""clones""=ordered_clones))1189}1190 1191#' @title get_evofreq_labels1192#' @param freq_frame Properly formatted dataframe from \code{\link{get_evofreq}}.1193#' @param evofreq_plot Plot returned from \code{\link{plot_evofreq}}.1194#' @param apply_labels Whether to plot the labels. Must include evofreq_plot if TRUE. Default FALSE.1195#' @param custom_label_text Custom labels to use (e.g. c(""Mutant X"",""Mutant Y""))1196#' @param extant_only Get labels for only the extant clones.1197#' @param clone_list Only label clones with this id. Passed as a vector.1198#' @param line_color Single color value. Default ""black"".1199#' @param adj.factor Value to scale the label and lines by so that they are far enough from the plot. Default value is 10 (a 10th of the distance between max/min at x from max/min across x).1200#' @return Dataframe containing the necessary coordinates to add labels in the proper locations.