Short term change #94

Daenarys8 · 2024-10-09T06:25:28Z

Function to return short term change in a dataframe

Signed-off-by: Daena Rys <[email protected]>

Daenarys8 · 2024-10-10T04:08:33Z

@antagomir @TuomasBorman

TuomasBorman · 2024-10-10T04:30:40Z

In holiday currently. I will check this early next week.

Signed-off-by: Daena Rys <[email protected]>

TuomasBorman

This PR is probably also not the most urgent ones

The calcultion part / core functionality is rather simple. The function should:

Get data with mia::meltSE (transformations sould be applied in prior to this function)
Calculate short term change based on group and time.col
Return the data.frame

After this, user could plot the results if wanted

R/shortTermChange.R

antagomir · 2024-11-03T20:24:26Z

Good to finalize this as proposed.

antagomir · 2024-11-03T20:29:18Z

Example shows the following:

# Load time series data
#' data(minimalgut)
#' tse <- minimalgut
#' 
#' short_time_labels <- c("74.5h", "173h", "438h", "434h", "390h")
#' 
#' # Subset samples by Time_lable and StudyIdentifier
#' tse <- tse[, !(colData(tse)$Time_label %in% short_time_labels)]
#' tse <- tse[, (colData(tse)$StudyIdentifier == "Bioreactor A")]
#' 
#' # Plot short term change in abundance
#' shortTermChange(tse, rarefy = TRUE, plot = TRUE) + ggtitle("Bioreactor A")

Is this mean to be an analysis or visualization function?

The example could explain a bit more what this aims to show. I am not sure if I understand the motivation of this function from the manpage description.

It would be good to have a separate function, like getShortTermChange / addShortTermChange and then a separate function for visualization (unless it is easy to do with existing functions, like miaViz::plotSeries).

Is this linked to some miaTime issue? Issue #7 ?

Signed-off-by: Daena Rys <[email protected]>

Daenarys8 · 2024-11-08T11:50:53Z

#'
#' # Plot short term change in abundance
#' shortTermChange(tse, rarefy = TRUE, plot = TRUE) + ggtitle("Bioreactor A")
Is this mean to be an analysis or visualization function?

original reference had plotting as part of the function. It is now removed and the function returns a dataframe with short term changes

Is this linked to some miaTime issue? Issue #7 ?

yes, syncomR consumed phyloseq object.

Signed-off-by: Daena Rys <[email protected]>

Daenarys8 · 2024-11-08T14:43:26Z

Moreover given the nature of the data frame results, I am not sure plotSeries is the right visualization method for it. Suppose we somehow addShortTermChange using one of the tse slots(say maybe reducedDims), it would still require it's own ggplot wrapper .
I just checked but to be sure @TuomasBorman do we have miaViz function that takes dataframes?
For reference, this is how the plots from syncomR could look like linear or polarized:

and of course we can also discuss how we want to represent our own plots.

TuomasBorman

Looks nice

R/getShortTermChange.R

TuomasBorman · 2024-11-08T17:45:48Z

Moreover given the nature of the data frame results, I am not sure plotSeries is the right visualization method for it. Suppose we somehow addShortTermChange using one of the tse slots(say maybe reducedDims), it would still require it's own ggplot wrapper . I just checked but to be sure @TuomasBorman do we have miaViz function that takes dataframes? For reference, this is how the plots from syncomR could look like linear or polarized:

and of course we can also discuss how we want to represent our own plots.

There is no that kind of restriction that miaViz cannot have functions that take df as input. In fact, loadings plotting function can take df as input,

This outputs currently data.frame and I think it makes sense. reducedDims is not correct because this does not reduce dimensionality. This could be stored to metadata and then we could have a plotting function for TreeSE and df.

antagomir · 2024-11-10T10:00:48Z

Or can it be part of colData? Could be simpler. This is a measure that is available on a per-sample basis. Except that some samples may have NAs if they are missing the comparison.

TuomasBorman · 2024-11-10T10:17:05Z

Or can it be part of colData? Could be simpler. This is a measure that is available on a per-sample basis. Except that some samples may have NAs if they are missing the comparison.

That could be possible but then columns would denote features. It might be confusing and there might be lots of features to add in the columns

antagomir · 2024-11-10T10:21:30Z

Or can it be part of colData? Could be simpler. This is a measure that is available on a per-sample basis. Except that some samples may have NAs if they are missing the comparison.

That could be possible but then columns would denote features. It might be confusing and there might be lots of features to add in the columns

Yep, not so optimal.

antagomir · 2024-11-10T10:32:29Z

"This function essentially calculates the change in abundance over time for a microbe by abundancet2/ abundancet1."

We already have getStepwiseDivergence?

getStepwiseDivergence calculates x(t + d) - x(t) for abundance x at time points t+d vs. t any user defined interval d along the time series.

getShortTermChange calculates this per taxa. One can generalize and use it for a more general interval x(t+d)/x(t), analogous to getStepwiseDivergence.

getShortTermChange suggests to do x(t+d)/x(t) instead of x(t+d)-x(t). But both are available and the choice may also depend on whether the data is logarithmic. log(a/b) = log(a)-log(b) although zeroes may cause trouble. It could be enough to implement the difference "x-y" and the rest ("x/y") can be always dealt (by the user?) with log/exp operations the data if necessary.

Signed-off-by: Daena Rys <[email protected]>

TuomasBorman

Looks good, some suggestions

TuomasBorman · 2024-11-14T11:00:03Z

R/getShortTermChange.R

+    # Reshape data and calculate growth metrics
+    assay_data <- meltSE(x, assay.type = assay.type, 
+                         add.col = time.col, row.name = "Feature_ID")
+    assay_data <- assay_data %>%


Use |> operators. They are used in elsewhere

TuomasBorman · 2024-11-14T11:00:26Z

R/getShortTermChange.R

+#' @export
+setMethod("addShortTermChange", signature = c(x = "SummarizedExperiment"),
+          function(x, assay.type = "counts", name = "short_term_change", ...){
+              # Calculate short term change
+              res <- getShortTermChange(x, ...)


Check BiocHCeck::BiocHeck() at least indentation is off

TuomasBorman · 2024-11-14T11:00:57Z

R/getShortTermChange.R

+#'          \url{https://github.com/nwisnoski/ul-seedbank}. Their approach is based on
+#'          the following calculation log(present abundance/past abundance).
+#'          Also a compositional version using relative abundance similar to
+#'          Brian WJi, Sheth R et al
+#'          \url{https://www.nature.com/articles/s41564-020-0685-1} can be used.
+#'          This approach is useful for identifying short term growth behaviors of taxa.


Remove this indentation

And also from above

TuomasBorman · 2024-11-14T11:02:03Z

R/getShortTermChange.R

+#' 
+#' short_time_labels <- c("74.5h", "173h", "438h", "434h", "390h")
+#' 
+#' # Subset samples by Time_label and StudyIdentifier


Instead of this comment, say what it is the goal: get specified time points from certain bioreactor

TuomasBorman · 2024-11-14T11:05:00Z

R/getShortTermChange.R

+#' \code{\link[SummarizedExperiment:SummarizedExperiment-class]{SummarizedExperiment}}
+#' object with these results in its \code{metadata}.
+#' 
+#' @details This approach is used by Wisnoski NI and colleagues


You could use math notation, it would be easier to read. Also the calculated values include other measurements than this

You could describe them also

TuomasBorman · 2024-11-14T11:23:58Z

tests/testthat/test-getShortTermChange.R

+    # Should still return a dataframe
+    short_time_labels <- c("74.5h", "173h", "438h", "434h", "390h")
+    # Subset samples by Time_label and StudyIdentifier 
+    tse_filtered <- tse[, !(tse$Time_label %in% short_time_labels)]
+    tse_filtered <- tse_filtered[, (tse_filtered$StudyIdentifier == "Bioreactor A")]
+
+    expect_true(all(!(tse_filtered$Time_label %in% short_time_labels)))
+
+    result <- getShortTermChange(tse_filtered, time.col = "Time.hr")
+    # Expected output is a dataframe
+    expect_true(is.data.frame(result))  
+    expect_true("growth_diff" %in% colnames(result))
+    # Test some expected properties (e.g., that growth_diff isn't all NAs)
+    expect_false(all(is.na(result$growth_diff)))


The values should be also checked

TuomasBorman · 2024-11-14T11:25:45Z

R/getShortTermChange.R

+setMethod("getShortTermChange", signature = c(x = "SummarizedExperiment"),
+    function(x, assay.type = "counts", ...){


I think the function should also take into account the grouping, i.e., from which patient/bioreactor the data is coming

TuomasBorman · 2024-11-14T12:05:42Z

R/getShortTermChange.R

+        arrange( !!sym(time.col) ) %>%
+        group_by(Feature_ID) %>%
+        mutate(
+            time_lag = !!sym(time.col) - lag( !!sym(time.col) ), 
+            growth_diff =!!sym(assay.type) - lag(!!sym(assay.type)),
+            growth_rate = (!!sym(assay.type) - lag(!!sym(assay.type))) / lag(!!sym(assay.type)),
+            var_abund = (!!sym(assay.type) - lag(!!sym(assay.type))) / time_lag
+        )
+    return(assay_data)


I think this should take into account duplicated values. I tink you could give warning first and then try something like this

df <- meltSE(tse, add.col = TRUE) df |> group_by(StudyIdentifier, FeatureID, Time.hr) |> mutate(mean_abundance_for_timepoint = mean(counts, na.rm = TRUE)) |> distinct(StudyIdentifier, FeatureID, Time.hr, .keep_all = TRUE) |> ungroup() |> arrange(StudyIdentifier, FeatureID, Time.hr) |> group_by(StudyIdentifier, FeatureID) |> mutate(growth_diff = mean_abundance_for_timepoint - lag(mean_abundance_for_timepoint)) |> ungroup()

TuomasBorman · 2024-11-14T12:07:30Z

R/getShortTermChange.R

+    grs.all$Feature_IDabb <- toupper(abbreviate(grs.all$Feature_ID, 
+                                               minlength = 3, 
+                                               method = "both.sides"))


These feature names come directly from rownames()? Then I think they should be kept unchanged

TuomasBorman · 2024-11-14T12:08:06Z

R/getShortTermChange.R

+        ########################### Growth Metrics ############################
+        grwt <- .calculate_growth_metrics(x, assay.type = assay.type, ...)
+        # Clean and format growth metrics
+        grs.all <- .clean_growth_metrics(grwt, ...)


Not sure if these are needed, as these are quite basic metric that can be easily calculated

Daenarys8 added 5 commits August 26, 2024 09:34

short term change

234bc47

Signed-off-by: Daena Rys <[email protected]>

up

4e9568a

Signed-off-by: Daena Rys <[email protected]>

add plot option

19937c7

Signed-off-by: Daena Rys <[email protected]>

update fxn

e451d86

Signed-off-by: Daena Rys <[email protected]>

tests for short term change

07e33b4

Signed-off-by: Daena Rys <[email protected]>

Daenarys8 added 2 commits October 10, 2024 11:46

import functions

452a709

Signed-off-by: Daena Rys <[email protected]>

update dependencies

6755ca9

Signed-off-by: Daena Rys <[email protected]>

TuomasBorman requested changes Oct 14, 2024

View reviewed changes

Daenarys8 added 2 commits November 8, 2024 09:57

Merge branch 'devel' into shortTermChange

7501031

update short term change

d0ac31f

Signed-off-by: Daena Rys <[email protected]>

update

99b98e3

Signed-off-by: Daena Rys <[email protected]>

TuomasBorman requested changes Nov 8, 2024

View reviewed changes

Daenarys8 added 3 commits November 13, 2024 12:47

update

9529271

Signed-off-by: Daena Rys <[email protected]>

update

b2afe4b

Signed-off-by: Daena Rys <[email protected]>

update

c34b0d7

Signed-off-by: Daena Rys <[email protected]>

Daenarys8 requested a review from TuomasBorman November 14, 2024 10:29

TuomasBorman requested changes Nov 14, 2024

View reviewed changes

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Short term change #94

Short term change #94

Daenarys8 commented Oct 9, 2024

Daenarys8 commented Oct 10, 2024

TuomasBorman commented Oct 10, 2024

TuomasBorman left a comment

antagomir commented Nov 3, 2024

antagomir commented Nov 3, 2024 •

edited

Loading

Daenarys8 commented Nov 8, 2024

Daenarys8 commented Nov 8, 2024

TuomasBorman left a comment

TuomasBorman commented Nov 8, 2024

antagomir commented Nov 10, 2024

TuomasBorman commented Nov 10, 2024

antagomir commented Nov 10, 2024

antagomir commented Nov 10, 2024

TuomasBorman left a comment

TuomasBorman Nov 14, 2024

TuomasBorman Nov 14, 2024

TuomasBorman Nov 14, 2024

TuomasBorman Nov 14, 2024

TuomasBorman Nov 14, 2024

TuomasBorman Nov 14, 2024

TuomasBorman Nov 14, 2024

TuomasBorman Nov 14, 2024

TuomasBorman Nov 14, 2024

TuomasBorman Nov 14, 2024

TuomasBorman Nov 14, 2024

TuomasBorman Nov 14, 2024

		setMethod("getShortTermChange", signature = c(x = "SummarizedExperiment"),
		function(x, assay.type = "counts", ...){

Short term change #94

Are you sure you want to change the base?

Short term change #94

Conversation

Daenarys8 commented Oct 9, 2024

Daenarys8 commented Oct 10, 2024

TuomasBorman commented Oct 10, 2024

TuomasBorman left a comment

Choose a reason for hiding this comment

antagomir commented Nov 3, 2024

antagomir commented Nov 3, 2024 • edited Loading

Daenarys8 commented Nov 8, 2024

Daenarys8 commented Nov 8, 2024

TuomasBorman left a comment

Choose a reason for hiding this comment

TuomasBorman commented Nov 8, 2024

antagomir commented Nov 10, 2024

TuomasBorman commented Nov 10, 2024

antagomir commented Nov 10, 2024

antagomir commented Nov 10, 2024

TuomasBorman left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

antagomir commented Nov 3, 2024 •

edited

Loading