R/summarise_lipids.R

In GitShepard1/lipidtools:

#' @title parsing the raw experimental data
#'
#' @description  This combines head with a number of columns
#' @param object takes an object generated by parse_raw function
#' @param grouping declare the grouping variable to summarise by e.g. CHAIN, SATURATION, REGION, etc
#' @param statistic choose summary statistic e.g. 'means +/- se', 'response' for difference between treatments, and 'percentage' for percentage out of total concentration
#' @param region_order declare region order
#' @import tidyverse
#' @import purrr
#' @import magrittr
#' @import broom
#' @export summarise_lipids
#' @return NULL


setClass("lipidSummary", representation(treatSummary = "tbl",
                                        groupSummary = 'tbl',
                                        grouping = "character",
                                        statistic = 'character',
                                        results = 'tbl',
                                        regions = 'character',
                                        colours = "vector",
                                        xvariable="character",
                                        xorder = "character",
                                        yvariable="character",
                                        ytitle="expression")
)

setMethod("show", "lipidSummary",
          function(object){

            cat('Summary Details \n')
            cat(" Statistic             :", object@statistic, "\n")
            cat(" Summarise by          :", object@grouping, '\n')
            cat(" Brain Regions         :", object@regions, '\n')


          })


summarise_lipids = function(object,
                            grouping,
                            statistic = 'means +/- se',
                            region_order = c("BLA","CeA","CB","FB","DH","VH")
){


  # parse_quosures for passing string arguments
  as_quosure <- function(strs) rlang::parse_quosures(paste(strs, collapse=";"))

  grouping_sel = unique(c(grouping, 'TREATMENT', 'REGION'))

  basedata = object@data %>%
    ungroup %>%
    group_by(!!!as_quosure(grouping_sel)) %>%
    dplyr::mutate(i = 1) %>%
    summarise(mean_sum = sum(FILT.MEAN, na.rm = T),
              se_sum = sum(FILT.SEM, na.rm = T),
              n = sum(i)
    ) %>%
    ungroup %>%
    sat_lookup(keys = grouping) %>%
    dplyr::mutate(treat = factor(paste(REGION, TREATMENT, sep ='\n')))

  # treatment summary
  treatsum = basedata %>%
    group_by(REGION, TREATMENT) %>%
    dplyr::summarise(
      mean_sum = sum(mean_sum, na.rm = T),
      se_sum = sum(se_sum, na.rm = T)
    ) %>%
    dplyr::mutate(treat = factor(paste(REGION, TREATMENT, sep ='\n'))) %>%
    arrange(match(REGION,  region_order),
            match(TREATMENT, object@treatments))



  # t tests

  if (length(object@treatments) == 2){

    statistic_df = object@data %>%
      dplyr::select(CHAIN, REGION, TREATMENT, paste0('M', 1:8,'')) %>%
      gather(SUBJECT_NO, VALUE, -CHAIN, -REGION, -TREATMENT) %>%
      group_by(REGION, TREATMENT, SUBJECT_NO) %>%
      dplyr::summarise(total = sum(VALUE)) %>%
      ungroup %>%
      group_by(REGION) %>%
      nest() %>%
      dplyr::mutate(ttest = map(data, function(df){

        t.test(total~TREATMENT, data = df) %>% tidy
      })) %>%
      unnest(ttest) %>%
      dplyr::mutate(stars = map(p.value, function(p) {


        if(p <=0.0001){

          "****"

        } else if(p <=0.001){

          "***"

        } else if(p <=0.01){

          "**"

        } else if(p <=0.05){

          "*"

        } else if(p >0.05){

          "ns"

        }
      })
      ) %>% unnest(stars) %>%
      dplyr::select(REGION, estimate, estimate1, estimate2, statistic, conf.low, conf.high, p.value, stars) %>%
      arrange(match(REGION,  region_order))

  } else {

    statistic_df = tibble()

  }


  # grouping summary
  groupsum = basedata %>%
    group_by(REGION, TREATMENT) %>%
    nest() %>%
    dplyr::mutate(total = map(data, function(df) sum(df$mean_sum))) %>%
    unnest(total) %>%
    unnest(data) %>%
    dplyr::mutate(percent = mean_sum/total)

  names(groupsum)[which(names(groupsum) == grouping)]<- 'FILL_BY'

  colour_pal = groupsum %>%
    distinct(FILL_BY, COLOUR, ORDER) %>%
    arrange(ORDER) %>%
    pull(COLOUR) %>%
    unique

  xaxis_order = groupsum %>%
    dplyr::select(TREATMENT, REGION) %>%
    arrange(match(REGION,  region_order),
            match(TREATMENT, object@treatments)) %>%
    dplyr::mutate(treat = paste(REGION, TREATMENT, sep ='\n')) %>%
    pull(treat) %>%
    unique

  # then
  if (statistic == 'response') {

    groupsum = groupsum %>%
      select(REGION, TREATMENT, FILL_BY, mean_sum, ORDER) %>%
      spread(TREATMENT, mean_sum) %>%
      dplyr::mutate(diff = get(object@treatments[1]) - get(object@treatments[2]))

  } else if (statistic == 'percentage') {

    groupsum = groupsum %>%
      dplyr::select(REGION, TREATMENT, FILL_BY, percent, treat, ORDER)


  } else if (statistic == 'means +/- se'){

    groupsum = groupsum %>%
      dplyr::select(REGION, TREATMENT, FILL_BY, mean_sum, se_sum, treat, ORDER)


  }


  return_class = new("lipidSummary",
                     treatSummary = treatsum,
                     groupSummary = groupsum,
                     grouping = grouping,
                     statistic = statistic,
                     results = statistic_df,
                     regions = region_order,
                     colours = colour_pal,
                     xvariable = ifelse(statistic == 'response', 'REGION', 'treat'),
                     xorder = if(statistic == 'response') region_order else xaxis_order,
                     yvariable = if(statistic == 'response') "diff" else if(statistic=='percentage') "percent" else "mean_sum",
                     ytitle = if(statistic != 'percentage') expression(paste(mu, "mol/mg tissue")) else expression("% total FFA")
  )

  return(return_class)

}