MSstats removes features with few measurements although removeFewMeasurements is set to False

[Romanos Siaperas]

Hi, I noticed that OpenMStoMSstatsFormat and dataProcess remove features with few measurements (in less than 3 runs) although the option is set to False. I have 3 replicates per condition so I want to keep features present in 2 runs. I am working with the output of quantms and MSstats v4.8.2

When reading the input file prior to MSstats I see features present in 1 or 2 runs.

> data <- read.csv("secretome/quantms/noMBR/proteomicslfq/sdrf_openms_design_msstats_in.csv", header = TRUE, sep = ',')

> data %>% group_by(PeptideSequence, PrecursorCharge) %>% summarise(RunCount = sum(!is.na(Intensity))) %>% ungroup() %>% count(RunCount) %>% mutate(Frequency = n / sum(n))

`summarise()` has grouped output by 'PeptideSequence'. You can override using the `.groups` argument.

# A tibble: 6 × 3

RunCount n Frequency 

<int> <int> <dbl>
1 1 2000 0.380

2 2 978 0.186

3 3 125 0.0237
4 4 374 0.0710
5 5 659 0.125

6 6 1128 0.214

Now if I read this in msstats

> data_msstats <- OpenMStoMSstatsFormat(data, useUniquePeptide = TRUE, removeFewMeasurements=FALSE, removeProtein_with1Feature = FALSE, + summaryforMultipleRows = max, use_log_file = FALSE)

INFO [2023-06-02 12:32:34] ** Raw data from OpenMS imported successfully.

INFO [2023-06-02 12:32:36] ** Raw data from OpenMS cleaned successfully.

INFO [2023-06-02 12:32:36] ** Using annotation extracted from quantification data. INFO [2023-06-02 12:32:36] ** Run labels were standardized to remove symbols such as '.' or '%'.

INFO [2023-06-02 12:32:36] ** The following options are used:

- Features will be defined by the columns: PeptideSequence, PrecursorCharge, FragmentIon, ProductCharge

- Shared peptides will be removed.

- Proteins with single feature will not be removed

. - Features with less than 3 measurements across runs will be kept.

INFO [2023-06-02 12:32:37] ** Features with all missing measurements across runs are removed.

INFO [2023-06-02 12:32:37] ** Shared peptides are removed.

INFO [2023-06-02 12:32:37] ** Multiple measurements in a feature and a run are summarized by summaryforMultipleRows: max

INFO [2023-06-02 12:32:37] ** Features with all missing measurements across runs are removed.

INFO [2023-06-02 12:32:37] ** Run annotation merged with quantification data. INFO [2023-06-02 12:32:37] ** Features with one or two measurements across runs are removed.

INFO [2023-06-02 12:32:37] ** Fractionation handled.

INFO [2023-06-02 12:32:38] ** Updated quantification data to make balanced design. Missing values are marked by NA INFO [2023-06-02 12:32:38] ** Finished preprocessing. The dataset is ready to be processed by the dataProcess function.

> data_msstats %>% group_by(PeptideSequence, PrecursorCharge) %>% summarise(RunCount = sum(!is.na(Intensity))) %>% ungroup() %>% count(RunCount) %>% mutate(Frequency = n / sum(n))

`summarise()` has grouped output by 'PeptideSequence'. You can override using the `.groups` argument.

# A tibble: 4 × 3

RunCount n Frequency

<int> <int> <dbl>

1 3 125 0.0547

2 4 374 0.164

3 5 659 0.288

4 6 1128 0.493

So MSstats only kept features present in >3 runs. The log output says first that

Features with less than 3 measurements across runs will be kept.

and then

INFO  [2023-06-02 12:32:37] ** Features with one or two measurements across runs are removed.

I tried the dataProcess function giving as input the data and not the data_msstats dataframe but the same filtering is applied again and prints that these features are removed.

> processed.quant <- MSstats::dataProcess(raw = data, + normalization = "equalizeMedians", + featureSubset = "all", + summaryMethod = 'TMP', + censoredInt = "NA", + MBimpute = TRUE, + use_log_file = FALSE)

INFO [2023-06-02 12:35:40] ** Features with one or two measurements across runs are removed.

INFO [2023-06-02 12:35:40] ** Fractionation handled.

INFO [2023-06-02 12:35:40] ** Updated quantification data to make balanced design. Missing values are marked by NA

INFO [2023-06-02 12:35:40] ** Log2 intensities under cutoff = 15.542 were considered as censored missing values.

INFO [2023-06-02 12:35:40] ** Log2 intensities = NA were considered as censored missing values.

INFO [2023-06-02 12:35:40] ** Use all features that the dataset originally has.

Checked again that there are no proteins with more than 3 missing values (I have 6 samples in total).

> sampleQuant <- quantification(processed.quant, use_log_file = FALSE)

> table(rowSums(is.na(sampleQuant[, 2:7])))  

0   1   2   3

416 112  69  28

[Mateusz Staniak]

Hi,

thanks for reporting the issue, do you have an example subset of data to reproduce it?

IIRC there are two or three spots in the workflow where we remove features with few observations and it is likely that only one is covered by the parameter. I will check what we can do about it.

Kind regards
Mateusz

[Romanos Siaperas]

Thank you for your response. I attach my dataset.

Best regards,

Romanos

[Mateusz Staniak]

Hi,

it seems that in your case 1-2 observations per feature only occur in proteins that have a single feature. They are removed based on the value of "removeProtein_with1Feature" parameter.

So in general case, changing value of this parameter would work.

However, there is an additional problem. We also remove features with few measurements while handling fractions and that is not regulated by the removeFewMeasurements parameters.

It was intended behavior, however it's clear that it's not the best idea here. I will update MSstats to change this, but figuring out the best way to handle it might takes a few more days.

You can either:

- wait for the update,

- run OpenMStoMSstatsFormat with MSstatsConvert extended workflow and change parameter "handle_fractions" in MSstatsBalancedDesign function to FALSE. This means pretty much running the internal code of OpenMStoMSstatsFormat with just that one parameter modified.

Kind regards,
Mateusz

[Romanos Siaperas]

Hi Mateusz,

I had already set the "removeProtein_with1Feature" parameter to False but these features are still being removed. I tried setting "handle_fractions" = FALSE in the OpenMStoMSstatsFormat but I still get the same filtered output as before, and it removes measurements present in a few runs.

The output of the quantms pipeline is ready for the dataProcess function, so I wouldn't mind if this filtering was only applied in the OpenMStoMSstatsFormat function. However, when running the dataProcess function using the output of quantms (which contains features present in 1 or 2 runs), it still removes features with few measurements.

Therefore, I will wait for the update.

Thank you very much for your time and effort.

Best regards,

Romanos

[Romanos Siaperas]

Dear Mateusz,

I wanted to follow up regarding the update. I understand that these things take time, but I was wondering if there is any rough estimate for when the update might be released. The reason I ask is that I need to schedule the downstream analyses for my project, and having a draft estimation of the waiting time would be quite helpful for planning purposes.

Thank you once again for your assistance.

Best regards,
Romanos

[Mateusz Staniak]

Hi,

unfortunately it requires a big change in MSstats workflow, so I didn't push any changes. Later this week I will provide a short piece of code that fixes this issue using MSstatsConvert workflow and I will consider making an update to MSstats, but it's non-trivial.

Kind regards,

Mateusz

[weixian...@gmail.com]

Hi Mateusz,

I read through this thread, I also ran into similar situation. 

my experimental design is that I have two conditions, each with duplicates, since it's immunopreicipation samples, negative controls are expected to have some completely missing proteins. There is one known interacting protein, which has two peptides/features quantified in both replicates of the bait IP condition but none in the control IP condition, therefore, each feature only got 2 measurements and got removed by the first step of dataProcess, removingfewmeasurement, and I didn't notice where to turn it off. 

Some convention also filter out things observed in less than 70% of runs, but in IP settings, let's say I have more than just bait vs. negative control IP, but several differentially treated bait IP conditions, while only one of the conditions having measurements for a feature, in a duplicate setting, it would be only two measurements, and by the rule of this convention and also MSStats dataProcess function, the feature would be removed, from biological point of view, it wouldn't make sense.

Combining both of our cases, to keep features having measurements in over 50% of the runs from the same condition seems make the most sense, or impute the rest of runs in this case. Wondering if there's any option that I wasn't aware of to fix my problem, or should I also need to expect the such fix in the future update?

Best,

Weixian 

[Mateusz Staniak]

Hi,

for OpenMS, please temporarily replace converter call by the following code:

(input is the OpenMS file)

input = MSstatsConvert::MSstatsImport(list(input = input),
                                      "MSstats", "OpenMS")
input = MSstatsConvert::MSstatsClean(input)
annotation = MSstatsConvert::MSstatsMakeAnnotation(input, NULL)

feature_columns = c("PeptideSequence", "PrecursorCharge",  "FragmentIon", "ProductCharge")
input = MSstatsConvert::MSstatsPreprocess(
  input,
  annotation,
  feature_columns,
  remove_shared_peptides = TRUE,
  remove_single_feature_proteins = F,
  feature_cleaning = list(
    remove_features_with_few_measurements = FALSE,
    summarize_multiple_psms = max))
colnames(input)
input[, Fraction := 1]

input = MSstatsConvert::MSstatsBalancedDesign(input, feature_columns, handle_fractions = FALSE)

Answering both questions, I will update MSstats to stop it from removing features in scenarios like this.

Kind regards,

Mateusz

[weixian...@gmail.com]

Hi Mateusz,

Thanks for posting the solution for OpenMS converter. Just to follow up with this issue, have you had any of PD/MQ/Philosopher converter updated to accommodate the afore-discussed cases? If I want to use quant tools other than OpenMS, could I adapt the code you wrote for OpenMS converter?

Best,

Weixian

[Mateusz Staniak]

Hi,

I will push a general fix to this issue to Bioconductor later this week. After this update, the parameter "removeFewMeasurements" will properly control this behavior.

Kind regards,

Mateusz Staniak

[Mateusz Staniak]

Hi,

I pushed this update to Bioconductor and Github. It's live on Github and should be up soon on Bioc. Both MSstatsConvert and MSstats require new installation. Please let me know if your problem was solved.

Kind regards,
Mateusz

[weixian...@gmail.com]

Thanks for the update, will test it soon and let you know.

Best regards

Weixian

[weixian...@gmail.com]

Hi Mateusz,

Just realized that fragpipe msstats input doesn't have to go through a converter, while using dataProcess function, it's still affected by remove few measurements action (Features with one or two measurements across runs are removed). It results in proteins with two peptides faithfully quantified in both replicates of one condition but none in the two replicates of control got filtered out. Do you have any suggestion or fix on this issue?

Best,

Weixian

[Devon Kohler]

Hi Weixian,

We actually added a FragPipe converter that does some filtering that was missing from the original output. Can you try the converter and see  if it fixes this issue?

Best,

Devon

[weixian...@gmail.com]

Oh, I actually tried that, reported "MSstatsFragPipeFiles" is not a defined class. Might be my fragpipe version? Haven't dug into that error, my MSStats is v4.9.8, btw.

[Mateusz Staniak]

Hi,

please also update MSstatsConvert, that should solve the problem.

Kind regards,

Mateusz

[weixian...@gmail.com]

Tested, all works as expected.

Thanks a lot.

Weixian