PertEval-scFM: Benchmarking Single-Cell Foundation Models for Perturbation Effect Prediction

Clarifications

Mean baseline (eq. 6): The mean baseline is computed as the mean expression profile across all unperturbed cells for each gene. We will revise the manuscript to clarify this computation more explicitly.

Error bars in Figure 2: The large error bars reflect variability across splits and model replicates. Increasing the number of splits and replicates could reduce this variability, and we note that the error bars seem large in the context of the of the scFMs results as they perform similarly, comparing the errors in GEARS with the scFMs for the first spectral splits we see that they do not overlap. Thus, the overlapping errors are a reflection of the similar performance of the scFMs rather than the number of splits.

MSE per perturbation (Figure 3): The results reported in Figure 3 are on a per-split basis. Specifically, we averaged the results across foundation models per split. This is apparent from the legend in which we indicate through colour coding what data point belongs to what split.

Conclusion

We thank the reviewer for their thoughtful feedback, which has highlighted several areas for improvement. We have increased the power of our dataset by including Norman two-gene perturbations and are addressing dataset limitations by incorporating the Replogle dataset and plan to expand our benchmark further over time. Thank you again for your constructive comments.

Correlative investigations of model performance

We thank the reviewer for pointing our potential improvements to Figure 3. To address these concerns, we have updated Figure 3 (see updated manuscript, Figure 3a and 3b), where we have included the Spearman correlation coefficient to quantify these observed relationships and we have added two-gene perturbations to increase the statistical power of the analysis. The results indicate that prediction performance does degrade as E-distance increases, but the relationship is weak ($\text{Spearman’s R} = 0.13, \text{p} = 1.22 × 10^{⁻⁸}$). This suggests that the model’s capability to predict perturbation effects depends not only on the magnitude of the perturbation, but also on its distribution. For example, we observed that some pronounced perturbations, such as CEBPA, are well-predicted despite having a high E-distance. This shows how we can use the E-distance to pinpoint interesting perturbations and failure modes beyond perturbation strength.

The reviewer notes that the figure “seems to say outliers are harder to predict”, but this interpretation is not entirely accurate. As shown in Figure 3a, the model performed well on certain perturbations, such as CEBPA, despite being one of the most pronounced perturbations in the dataset, as measured by E-distance. Figure 4 provides insight into this apparent discrepancy, suggesting that the narrow peak of the effect inflates the E-distance while the long tail of minimal effects makes the overall distribution easier to predict.

Dataset splitting and SPECTRA

We appreciate the concern regarding the impact of dataset splitting. However, we believe SPECTRA’s contribution extends beyond what random splits or current experimental designs offer. In real-world applications like Perturb-seq, unpredictable biological distribution shifts are common. Random splits are not inherently easier or harder; a random split could align with a sparsification probability of 0.1 or 0.7, representing significantly different train-test dissimilarities. However, systematically quantifying this dissimilarity makes it possible to contextualise performance results and compare model robustness across different tasks or datasets. Random splits may therefore fail to capture these distribution shifts, often leading to overly optimistic assessments of model performance. SPECTRA addresses this by generating splits of increasing difficulty, starting with random splits and progressing to more challenging ones based on train-test similarity. This approach ensures that SPECTRA encompasses random splits but also evaluates a model’s robustness across a broader spectrum of biological variability. SPECTRA does not artificially inflate task difficulty but rather controls for the similarity between train and test split. This indeed increases test difficulty, but not artificially. While controlling for train-test similarity may lead to smaller datasets in harder splits, this is a reflection of the redundancy inherent to the dataset, not a limitation of SPECTRA. Moreover, SPECTRA is not just an evaluation tool but a valuable guide for experimental design. For example, in Perturb-seq, SPECTRA can reveal if a model generalises well within pathways but struggles across pathways, guiding researchers on which genes to test next. Random splits, by contrast, obscure such insights due to high train-test similarity. We argue that SPECTRA’s ability to highlight weaknesses in generalisability, simulate real-world distribution shifts, and inform experimental design makes it a high-impact addition to evaluation pipelines, especially for tasks like Perturb-seq where navigating distribution shifts is critical.

Biological relevance of loss-of-function simulations

We thank the reviewer for highlighting the potential limitations of simulating perturbations by nullifying the expression of the perturbed gene to generate scFM embeddings. While we acknowledge that this approach may appear counterintuitive from a biological perspective, it was chosen to address the challenges inherent in standardising perturbation simulations across foundation models. scFMs differ in gene expression representation. For instance, scGPT encodes gene expression as numerical values, whereas Geneformer uses rank-order representations. If we were to simulate upregulation by artificially increasing a gene's expression value in numerically encoded models, there is no biologically grounded way to determine the appropriate magnitude. Arbitrary modifications risk introducing biases that are inconsistent across models. For non-parametric models that use rank-order representations, simulating upregulation becomes even more ambiguous, as it is unclear how to re-rank genes to reflect such a perturbation.

On the other hand, nullifying the expression of a perturbed gene provides a uniform and unbiased approach that avoids model-specific inconsistencies, enabling fairer comparisons across diverse foundation models. While this approach is better suited to loss-of-function scenarios, we recognise its limitations in representing gain-of-function perturbations, such as those present in the Norman dataset. This limitation reflects a broader challenge in foundation model methodologies, which we briefly address in the discussion section of the paper: “key questions, such as how to represent perturbations in silico, $...$ need to be addressed” (line 469). Addressing these representational issues will be crucial for future development of foundation models to better capture the effects of both loss- and gain-of-function perturbations. Despite this, we believe that it is beyond the scope of our benchmarking paper to suggest alternative approaches for this issue. We hope our work will motivate further efforts to refine how perturbations are simulated to improve their biological relevance and predictive accuracy.

From the results presented in Tables 3 and 4, we observe a slight increase in prediction performance across certain models. However, the relative performance compared to baseline models diminishes as task complexity increases. This aligns with our prior observations that generalisability decreases for perturbations that are strong or atypical, which tend to be overrepresented in dual-gene perturbations. While the additional results provide further evidence of the trends observed in our initial analysis, they also underline the inherent limitations of current scFM embeddings in addressing the complexities of dual-gene perturbations and emphasise the necessity of future efforts to improve their robustness and biological relevance. Furthermore, we are currently working on including the Replogle et al. dataset in our benchmark. While compute and memory constraints currently pose challenges to processing this dataset, we are working on engineering solutions to implement this analysis and aim to report these results as soon as they are available.

We thank the reviewer for their detailed feedback and for recognising the strengths of our work, including the importance of the research question, the overall quality of the manuscript, the novelty of SPECTRA and contextual alignment, and our commitment to publishing negative findings. Below, we address the concerns and questions raised in detail.

Benchmarking on a single dataset We agree that the use of more datasets would broaden the scope of our analysis. We want to clarify that we originally focused on the Norman dataset because it has been shown to contain the strongest perturbation effect signal for genetic perturbations as measured by the E-distance [See Figure 3.b. in [1]] and is therefore best suited to our analysis. We hypothesised that if scFM embeddings contained biologically relevant information for perturbation effect prediction, it would be more likely to be contained in a dataset with a strong perturbation effect signal. Positive results in this context would indicate that biologically relevant perturbation effects can be predicted accurately with zero-shot scFM embeddings, while negative results would provide greater confidence that the observed limitations reflect true negatives, rather than being due to low-quality signal inherent to the dataset tested. In line with this argument, the Replogle dataset, despite containing a greater number of perturbations (~4,000 K562; ~200 RPE1), demonstrates much diminished perturbation effect signal overall [See Figure 3.b. in [1]]. We have expanded our analysis to include all the combinatorial perturbations in the Norman dataset, where originally we only considered single-gene perturbations. We present the updated results in Table 3 for the top 2000 HVGs and in Table 4 for the top 20 DEGs.

Table 3: Perturbation effect prediction evaluation across 2000 HVGs.

Table 4: Perturbation effect prediction evaluation across the top 20 DEGs.

Note: For two-gene perturbations split 0.5, there were not enough perturbations that passed our quality control to properly define the split.

[1] Peidli et al., 2024. scPerturb: harmonized single-cell perturbation data. https://www.nature.com/articles/s41592-023-02144-y

To implement the double-gene perturbations, we calculated the co-expression matrix similar to how we did for the single-gene perturbation, and averaged to obtain a unified co-expression matrix.

We have now updated the manuscript to reflect this, along with the following changes:

• Added dual-gene perturbations into Table 2, for evaluation across 2,000 Highly Variable Genes (HVGs).
• Added dual-gene perturbations into Table 3, for evaluation across 20 Differentially Expressed Genes (DEGs).
• Updated Figure 3 to include Spearman correlation and added all dual-gene perturbations to increase the power and statistical significance of the observed trend.
• Included explanation on dual-gene in silico perturbation representation in Section 2.1.2 and Appendix C1.1 and C2.
• Included interpretation of the dual-gene results in Sections 3.1, 3.2 and 3.3.
• Mentioned double-gene perturbation evaluation in the conclusion.

Changes have been marked in red, so they should be easy to spot.

is there only a single train-test split per sparsification probability?

There are three train-test splits per sparsification probability ($s$) because we run the experiments in triplicate with different seeds. To determine these splits, we construct a similarity graph where similar perturbations are connected. $s$ then controls how many edges are stochastically removed from the similarity graph to obtain subgraphs from which train and test perturbations can be sampled. This enforces a certain degree of cross-split overlap between train and test sets such that as $s$ increases, the train and test sets become more dissimilar. We explicitly quantify this cross-split overlap and illustrate its decline as $s$ increases (Appendix Figure E1.b). Figure 3 shows the splits across all replicates. Further details on SPECTRA, our train-test splitting procedure, and the role of the sparsification probability are provided in Section 2.2.3 and Appendix E of the manuscript.

And every perturbed gene only appears in a single split and sparsification probability?

A single perturbed gene can theoretically appear across multiple splits. However, SPECTRA ensures that the same or sufficiently similar perturbations do not simultaneously appear in both the train and test sets of any split. This separation minimises information leakage between the train and test sets and allows us to simulate distribution shifts.

how did you pick which sparsification probability to assign to a gene in Figure 3?

SPECTRA tracks the assignment of each perturbation to its respective train-test split. We identify a split by its sparsification probability, which runs from 0.0 to 0.7. For Figure 3 in particular, we are visualising the test perturbations for all splits and from SPECTRA we can link exactly which perturbation belongs to which train-test split. This can be appreciated in Figure 3a through the colour coding, where each test perturbation is assigned to a particular split, and in Figure 3b, we have organized all the test perturbations per split.

Thank you for bringing up these points. We hope these clarifications address your questions.

dual gene results

Thanks for adding this! Could you provide details on how you implemented the prediction, ideally in an updated version of the manuscript? Thanks!

MSE per perturbation (Figure 3): The results reported in Figure 3 are on a per-split basis. Specifically, we averaged the results across foundation models per split. This is apparent from the legend in which we indicate through colour coding what data point belongs to what split.

I'm still not sure I understand this: is there only a single train-test split per sparsification probability? And every perturbed gene only appears in a single split and sparsification probability? Else, how did you pick which sparsification probability to assign to a gene in Figure 3?

We sincerely thank the reviewer for acknowledging our extensive efforts to include two-gene perturbations across all models and baselines, including the newly added GEARS, evaluated under both HVG and DEG setups.

We appreciate the reviewer’s continued concern regarding the use of a single dataset. To address this, we are actively working on integrating the Replogle et al. (2022) dataset $1$ . Our preprocessing and featurisation pipeline is fully operational, and model training is currently underway. We anticipate completing the full evaluation within approximately one week. As an update, we provide preliminary results for three scFMs and two baselines across five SPECTRA splits, which highlight the progress we have made in broadening our dataset scope.

Regarding the SPECTRA splitting strategy, we respectfully maintain that its value extends beyond merely increasing task difficulty. SPECTRA provides a systematic and fair approach to evaluating the generalisability of methods under varying levels of distribution shift, a critical aspect in the context of single-cell datasets, which are particularly prone to batch effects and noisy measurements. The importance of assessing robustness to such shifts is further highlighted by the performance of task-specific models like GEARS. While GEARS outperforms other methods in general, it still exhibits a clear degradation in performance as sparsification probability increases. Furthermore, differences in the performance of embeddings from individual scFM models further underscore the relevance of this analysis. For example, scGPT experiences catastrophic performance degradation at a sparsification probability of 0.7, whereas UCE maintains performance levels comparable to the baseline (Appendix H.3). These insights are valuable for understanding the robustness and limitations of the evaluated models and are independent of the inherent complexity of the task itself. We hope this perspective clarifies the importance of SPECTRA in the broader context of model evaluation.

Thank you once again for your thoughtful feedback.

$1$ Replogle, J. M. et al. Mapping information-rich genotype–phenotype landscapes with genome-scale Perturb-seq. Cell 185, 2559–2575 (2022).

We thank the reviewer for their thoughtful feedback and for recognizing the strengths of our work, including the creation of a standardised evaluation framework, the inclusion of complementary and novel metrics, the evaluation of major scFMs with appropriate baselines, and the constructive paths forward suggested by our distribution shift analysis. Below, we address the reviewer’s concerns in detail:

Limited dataset scope

Regarding the use of additional datasets, we agree this would broaden the scope of our analysis. We want to clarify that we originally focused on the Norman dataset because it has been shown to contain a higher-quality and stronger perturbation signal $See Figure 3.b. in [1](https://www.nature.com/articles/s41592-023-02144-y/figures/3)$ and was therefore better suited to our analysis. We hypothesised that if scFM embeddings contained biologically relevant information for perturbation effect prediction, it would be more likely to manifest in a dataset with a strong signal like Norman. Positive results in this context would indicate that biologically relevant perturbation effects can be predicted accurately with zero-shot scFM embeddings, while negative results would provide greater confidence that the observed limitations reflect true negatives, rather than being due to low-quality signal inherent to the dataset tested.

In line with this argument, the Replogle dataset, despite containing a greater number of perturbations (~4,000 K562; ~200 RPE1), demonstrates a lower-quality perturbation signal overall $[1](https://www.nature.com/articles/s41592-023-02144-y)$ . This may be due in part to the complexity of its perturbation space: seeing more perturbations simultaneously may introduce signal that is diffuse and not conclusively tied to any single perturbation, resulting in increased noise. This characteristic makes it less well-suited for drawing strong conclusions about scFM embeddings in their current form. However, we agree with the reviewers that incorporating additional datasets like Replogle would strengthen the paper and increase the credibility of our conclusions. While compute and memory constraints currently pose challenges to processing the Replogle dataset, we are working on engineering solutions to enable this analysis. We are committed to reporting these results as soon as they are available.

[1] Peidli et al., 2024. scPerturb: harmonized single-cell perturbation data. https://www.nature.com/articles/s41592-023-02144-y

Methodological choices around MLP probes

We selected simple MLP probes with one hidden layer to isolate and evaluate the intrinsic information content of scFM embeddings. This minimalistic architecture reduces confounding effects from complex prediction heads and focuses on the utility of the embeddings themselves.

Sensitivity to architecture choices: Our experiments varying the number of layers and hidden dimensions (see Table 1) showed no significant improvement in performance. This suggests that the limitations lie in the embeddings themselves rather than the probe architecture.

Minor Note on Table 1

This has now been changed to “Transformer” in Table 1, thank you for spotting that!

Conclusion

We appreciate the reviewer’s insightful comments, which have highlighted important areas for improvement. We are addressing the dataset scope concern by incorporating the Replogle dataset into our benchmark and expanding the GitHub resource over time. We have investigated the impact of MLP parameter disparities and found limited sensitivity to architecture choices. Thank you again for your constructive feedback.

Parameter count disparities in MLP probes

We appreciate the reviewer’s concern about the disparities in MLP parameter counts due to varying embedding dimensions across scFMs. We acknowledge that this could influence training dynamics and make fair comparison challenging. To alleviate this concern, we trained a range of MLPs with increasing parameter count on our raw expression data to verify the effect of parameter number on results. We report our findings in Table 1, where it can be seen the increase in parameters has no effect on the MSEs obtained. We also include additional results using scBERT and scFoundation embeddings as input, with increasing parameter count. This further strengthens our conclusion that, despite an increase in the expressiveness of the underlying MLP model, the scFM embeddings do not contain sufficient biologically relevant information from the pre-training stage to provide increased performance for this task.

Table 1: Train and Test set results with MLPs of increasing parameter count

We agree that normalising model capacity — for instance, by using fixed-capacity projection layers to map all embeddings to a standardised dimensionality before passing them into the MLP — could theoretically provide a more equitable comparison. However, the results in Table 1 suggest that such normalisation is unlikely to change the overall conclusions. Additionally, since we are primarily evaluating the intrinsic information content of the embeddings, projecting them to a different dimensionality might introduce confounding factors, potentially obscuring the original information encoded in the embeddings.

Despite the disparities in parameter count, our current approach maintains alignment with the standard practice of tailoring MLP architectures to the embedding size to ensure sufficient expressiveness for each model. This guarantees that no embedding is artificially constrained or unfairly advantaged due to mismatched architecture capacity.

I thank the reviewers for their detailed response and for adding additional information to the paper. Since my main problem has now been sufficiently addressed I have decided to raise my score.

We sincerely thank the reviewer for raising our score and for the constructive feedback.

To address the original concerns about dataset scope, we are actively working on incorporating the Replogle dataset [1] into our analysis. This dataset introduces approximately 2,000 perturbations, significantly expanding the scope of our evaluation. Our preprocessing and featurization pipelines for this dataset are fully operational. However, due to its size, generating the necessary features and training all models across all splits requires substantial computational resources. We anticipate completing this process within one week.

In the meantime, we are pleased to share preliminary results in the updated PDF for three foundation models and two baselines across five SPECTRA splits. These results demonstrate our progress and reflect our commitment to broadening the study’s dataset scope as originally suggested.

Thank you again for your valuable input and insightful comments.

[1] Replogle, J. M. et al. Mapping information-rich genotype–phenotype landscapes with genome-scale Perturb-seq. Cell 185, 2559–2575 (2022).

Additional baseline

In Table 1 we include results where we have implemented GEARS on the raw expression, using SPECTRA splits. GEARS obtains strong predictive performance, outperforming baseline models, consistent with its purpose built architecture specifically designed for Perturbation effect prediction. Similarly to other models, GEARS achieves better performance in prediction of overall perturbation effect compared to effect on the top 20 differentially expressed genes (DEGs). This points to an inherent limitation in the ability to capture fine-grained perturbation-specific responses (such as in highly regulated DEGs). We again observe an increase in MSE with increasing sparsification probability, a pattern consistent with PertEval findings. This shows that increasing train-test dissimilarity reduces model robustness, emphasising the difficulty of generalising to out-of-distribution perturbations. This underscores the need for robust evaluation frameworks that take distribution shift into account and highlights both current limitations and future directions for research. In tandem with the finding scFM embeddings do not provide meaningful information for perturbation prediction, these results advocate for the development of specialised models for perturbation prediction tasks, as opposed to adapting general-purpose scFMs, which is one of the conclusions of our paper.

Table 1: MSE for perturbation effect prediction with GEARS trained on raw expression data across increasing sparsification probabilities

Overall MSE

Top 20 DEGs

Novelty

While we acknowledge that some of the evaluation metrics and models included in our benchmark have been previously explored, PertEval-scFM introduces several key contributions that add value to the field:

• Systematic distribution shift evaluation: The inclusion of the SPECTRA framework allows us to evaluate scFM robustness under varying train-test dissimilarity, a critical aspect for real-world applications that involve significant domain shifts.
• Comprehensive metrics: By integrating AUSPC, E-distance, and contextual alignment, PertEval-scFM provides a more holistic view of model performance, focusing on robustness and sensitivity to distribution shifts rather than accuracy alone.
• Open-source framework: Our codebase is designed to be extendable, enabling the community to incorporate new datasets, models, and metrics, fostering ongoing development in this area.

These contributions aim to establish a standardised and reproducible evaluation framework for scFMs, addressing a gap in the field despite the previously studied task and models.

Setting / modelling paradigm

Thank you for the feedback on this aspect of our work. We would first like to note that even in the absence of fine-tuning, if a signal relative to perturbation effect prediction were present in the sc-embeddings, it should still be measurable in a zero-shot setting. If no performance improvement comes from foundation model embeddings and all observed performance gains are due to fine-tuning, this contradicts the guiding principle that these models learn useful biological representations during pre-training. MLP probes are commonly used as an approach to answer the question of information content of embeddings in NLP and CV, as they evaluate representation quality while removing confounding effects pertaining to task-specific prediction heads, which introduce inductive biases [2,3,4]. We also note that previous work has investigated the performance of fine-tuned versions of scGPT and scFoundation [5]. This study found that simple linear baselines still outperformed these models, suggesting that fine-tuning does not fully address their limitations. Because fine-tuning performance is addressed in other studies and because it fundamentally goes against our approach of establishing existing information content, we do not include fine-tuning in our study design. The findings we present therefore highlight that scFM embeddings have not learnt useful biological patterns pertaining to the perturbation prediction task. This in and of itself is an important finding.

[2] Radford et al., 2021. Learning Transferable Visual Models From Natural Language Supervision. https://arxiv.org/abs/2103.00020

[3] Tenney et al., 2024. What do you learn from context? Probing for sentence structure in contextualized word representations https://openreview.net/forum?id=SJzSgnRcKX

[4] Boiarsky et al., 2023. A Deep Dive into Single-Cell RNA Sequencing Foundation Models. https://www.biorxiv.org/content/10.1101/2023.10.19.563100v1

[5] Ahlmann-Eltze et al., 2024. Deep learning-based predictions of gene perturbation effects do not yet outperform simple linear methods. https://www.biorxiv.org/content/10.1101/2024.09.16.613342v4

We thank the reviewer for their thoughtful feedback and detailed comments on our work. We appreciate the recognition of the relevance of genetic perturbation prediction, the inclusion of multiple scFMs, the clarity of our manuscript and figures, and the open-source nature of our benchmarking framework. Below, we address the reviewer’s concerns in detail:

Datasets

Regarding the use of additional datasets, we agree this would broaden the scope of our analysis. We want to clarify that we originally focused on the Norman dataset because it has been shown to contain a higher-quality and stronger perturbation signal [See Figure 3.b. in [1]] and was therefore better suited to our analysis. We hypothesised that if scFM embeddings contained biologically relevant information for perturbation effect prediction, it would be more likely to manifest in a dataset with a strong signal like Norman. Positive results in this context would indicate that biologically relevant perturbation effects can be predicted accurately with zero-shot scFM embeddings, while negative results would provide greater confidence that the observed limitations reflect true negatives, rather than being due to low-quality signal inherent to the dataset tested.

In line with this argument, the Replogle dataset, despite containing a greater number of perturbations (~4,000 K562; ~200 RPE1), demonstrates a lower-quality perturbation signal overall [1]. This may be due in part to the complexity of its perturbation space: seeing more perturbations simultaneously may introduce signal that is diffuse and not conclusively tied to any single perturbation, resulting in increased noise. This characteristic makes it less well-suited for drawing strong conclusions about scFM embeddings in their current form. However, we agree with the reviewers that incorporating additional datasets like Replogle would strengthen the paper and increase the credibility of our conclusions. While compute and memory constraints currently pose challenges to processing the Replogle dataset, we are working on engineering solutions to enable this analysis. We are committed to reporting these results as soon as they are available.

The reviewer raises an important point about the potential for our embedding generation approach (line 133) to influence the quality of embeddings. Our approach involves generating perturbed cell embeddings by setting the expression of the perturbed gene to zero in all cells where it is expressed, effectively simulating a perturbation in silico. This methodology is consistent with approaches used in prior scFM papers. We acknowledge that this approach might not fully capture the complexity of biological perturbations, but it provides a consistent framework for benchmarking scFMs. Future work could explore alternative embedding generation strategies tailored to specific biological contexts.

[1] Peidli et al., 2024. scPerturb: harmonized single-cell perturbation data. https://www.nature.com/articles/s41592-023-02144-y

We thank the reviewer for their constructive feedback, which has highlighted important areas for improvement. We believe that our focus on zero-shot evaluation, systematic distribution shift analysis, and comprehensive metrics provides significant contributions to the field, while our ongoing efforts to include additional datasets and advanced baselines will address the limitations raised. By continuing to expand the PertEval-scFM framework and its accompanying GitHub resource, we provide a valuable and evolving benchmarking tool for the community. Thank you again for your detailed comments.

We appreciate the reviewer’s thoughtful feedback and the opportunity to improve our work. PertEval-scFM provides significant contributions in robust evaluation under distribution shifts, comprehensive metrics, and transparent benchmarking, all aimed at advancing scFM research. We plan to expand our framework to include additional datasets like Replogle, provide more biological insights into atypical perturbations, and explore the integration of scFM embeddings into advanced predictive models. These updates, along with continuous expansion of our accompanying GitHub resource, will ensure that PertEval-scFM remains a valuable and evolving resource for the community. Thank you again for your constructive feedback.

Absence of fine-tuning component

We understand the reviewer's concern regarding the fact that we are not fine-tuning the scFMs. We would first like to note that even in the absence of fine-tuning, if a signal relative to perturbation effect prediction were present in the sc-embeddings, it should still be measurable in a zero-shot setting. If no performance improvement comes from foundation model embeddings and all observed performance gains are due to fine-tuning, this contradicts the guiding principle that these models learn useful biological representations during pre-training. MLP probes are commonly used as an approach to answer the question of information content of embeddings in NLP and CV, as they evaluate representation quality while removing confounding effects pertaining to task-specific prediction heads, which introduce inductive biases [3, 5, 6]. We also note that previous work has investigated the performance of fine-tuned versions of scGPT and scFoundation [7]. This study found that simple linear baselines still outperformed these models, suggesting that fine-tuning does not fully address their limitations. Because fine-tuning performance is addressed in other studies and because it fundamentally goes against our approach of establishing existing information content, we do not include fine-tuning in our study design. The findings we present therefore highlight that scFM embeddings have not learnt useful biological patterns pertaining to the perturbation prediction task. This in and of itself is an important finding.

[3] Radford et al., 2021. Learning Transferable Visual Models From Natural Language Supervision. https://arxiv.org/abs/2103.00020

[5] Tenney et al., 2024. What do you learn from context? Probing for sentence structure in contextualized word representations https://openreview.net/forum?id=SJzSgnRcKX

[6] Boiarsky et al., 2023. A Deep Dive into Single-Cell RNA Sequencing Foundation Models. https://www.biorxiv.org/content/10.1101/2023.10.19.563100v1

[7] Ahlmann-Eltze et al., 2024. Deep learning-based predictions of gene perturbation effects do not yet outperform simple linear methods. https://www.biorxiv.org/content/10.1101/2024.09.16.613342v4

scFM embeddings in GEARS and CELLOT

With regards to including results on scFM embeddings with GEARS and CellOT, we note that for GEARS we can only integrate gene embeddings per design, as GEARS explicitly uses a gene coexpression knowledge graph and Gene Ontology (GO)-derived relationships to learn embeddings that represent gene-specific perturbation effects. Cell embeddings, which encapsulate whole-cell states, would not align directly with these encoders unless significant architectural changes were made. This would involve revising the core design of GEARS to account for multiscale modelling, which we believe is beyond the scope of our paper. CellOT is designed to learn mappings between unpaired distributions of untreated and treated cells via Optimal Transport. In contrast, embeddings generated by scFM should capture intrinsic properties of cells without specific perturbation context. Integrating embeddings from zero-shot single-cell foundation models as inputs to CellOT will likely not align with CellOT's design, which requires raw single-cell measurements from both untreated and perturbed states to learn accurate mappings and embeddings. Therefore we believe it’s not appropriate for us to incorporate these prediction models within our benchmark.

MLP structure selection

We acknowledge the reviewer’s concern regarding the single-layer MLP probe and its potential limitations. However, our choice of MLP architecture was motivated by its simplicity, which allows for isolating the information content of scFM embeddings without introducing additional confounding factors from complex prediction heads. Indeed, one of the main aims of our work was to evaluate the representation quality of the scFM embeddings while minimising the amount of inductive bias introduced by model architecture. A simple, one-layer MLP structure has been used in other comparable studies [2,3]. This approach also allows us to evaluate representation robustness under distribution shift, without confounding results via task-specific prediction head inductive bias [4].

However, we also conducted experiments with varying MLP structures and found that increasing the complexity of the MLP had no significant impact on the results. We trained a range of MLPs with increasing parameter count on our raw expression data to verify the effect of parameter number on results. We report our findings in Table 1, where it can be seen the increase in parameters has no effect on the MSEs obtained. We also include additional results using scBERT and scFoundation embeddings as input, with increasing parameter count. This further strengthens our conclusion that, despite an increase in the expressiveness of the underlying MLP model, the scFM embeddings do not contain sufficient biologically relevant information from the pre-training stage to provide increased performance for this task.

Table 1: Train and Test set results with MLPs of increasing parameter count

Despite the disparities in parameter count, our current approach maintains alignment with the standard practice of tailoring MLP architectures to the embedding size to ensure sufficient expressiveness for each model. This guarantees that no embedding is artificially constrained or unfairly advantaged due to mismatched architecture capacity.

[2] Goyal et al., 2023. Finetune Like You Pretrain: Improved Finetuning of Zero-Shot Vision Models https://openaccess.thecvf.com/content/CVPR2023/html/Goyal_Finetune_Like_You_Pretrain_Improved_Finetuning_of_Zero-Shot_Vision_Models_CVPR_2023_paper.html

[3] Radford et al., 2021 Learning Transferable Visual Models From Natural Language Supervision https://arxiv.org/abs/2103.00020

[4] Heim et al. 2023, Towards Better Understanding of Domain Shift on Linear-Probed Visual Foundation Models https://openreview.net/forum?id=l188N6IZNY

We thank the reviewer for their feedback and constructive comments, as well as for acknowledging the strengths of our work, including the breadth of scFMs benchmarked, the inclusion of distribution shift scenarios, the comprehensive evaluation metrics, and the transparency and reproducibility of PertEval-scFM. Below, we address the reviewer’s concerns in detail:

Novelty and significance

We acknowledge that benchmarking scFMs on perturbation prediction tasks is not entirely novel. However, our work introduces unique contributions that differentiate PertEval-scFM: Distribution shift evaluation: we introduce systematic sparsification of train-test overlaps through the SPECTRA framework, providing a quantitative measure of how scFM performance is affected by increasing train-test dissimilarity. This approach offers insights into the robustness of scFM embeddings under realistic out-of-distribution scenarios, which are highly relevant to biological applications. Comprehensive metrics: the combination of AUSPC, E-distance, and contextual alignment provides a more detailed analysis of model performance than previous benchmarks. While existing studies have focused on accuracy-based metrics, we focus on robustness and sensitivity to distribution shifts, which are critical for real-world applications of scFMs. Regarding the use of additional datasets such as those in Wu et al. (2024), we agree this would broaden the scope of our analysis. We want to clarify that we originally focused on the Norman dataset because it has been shown to contain a higher-quality and stronger perturbation signal [See Figure 3.b. in [1]] and was therefore better suited to our analysis. We hypothesised that if scFM embeddings contained biologically relevant information for perturbation effect prediction, it would be more likely to manifest in a dataset with a strong signal like Norman. Positive results in this context would indicate that biologically relevant perturbation effects can be predicted accurately with zero-shot scFM embeddings, while negative results would provide greater confidence that the observed limitations reflect true negatives, rather than being due to low-quality signal inherent to the dataset tested.

In line with this argument, the Replogle dataset, despite containing a greater number of perturbations (~4,000 K562; ~200 RPE1), demonstrates a lower-quality perturbation signal overall [1]. This characteristic makes it less well-suited for drawing strong conclusions about scFM embeddings in their current form. However, we agree with the reviewers that incorporating additional datasets like Replogle would strengthen the paper and increase the credibility of our conclusions. While compute and memory constraints currently pose challenges to processing the Replogle dataset, we are working on engineering solutions to enable this analysis. We aim to report these results as soon as they are available. We appreciate the reviewer’s suggestion to provide a more quantitative and systematic analysis of strong/atypical perturbations and their biological context. To address this, we performed additional correlation analyses to investigate the relationship between the perturbation effect size, as measured by energy distance (E-distance), and model performance (MSE). In the updated Figure 3 (see the attached pdf), we have included the Spearman correlation coefficient to quantify these observed relationships and added two-gene perturbations to increase the statistical power of the analysis. The results indicate that prediction performance does degrade as E-distance increases, but the relationship is weak ($\text{Spearman’s R} = 0.13, \text{p} = 1.22 × 10^{⁻⁸}$). This suggests that the model’s capability to predict perturbation effects depends not only on the magnitude of the perturbation, but also on its distribution. For example, we observed that some pronounced perturbations, such as CEBPA, are well-predicted despite having a high E-distance. This can be explained by the narrow, pronounced peak of the effect (which inflates the E-distance) combined with a long tail of minimal effects, which makes the overall distribution easier to predict. This shows how we can use the E-distance to pinpoint interesting perturbations and failure modes beyond perturbation strength. This analysis reflects the role of PertEval-scFM as a benchmarking tool: the metrics we include (e.g. E-distance and AUSPC) are designed to enable model meaningful evaluation from a methodological as well as a biological perspective. They complement standard bioinformatic analyses approaches, by further connecting model predictions to meaningful biological mechanisms.

[1] https://www.nature.com/articles/s41592-023-02144-y

As stated in my previous comments, I expected to see more datasets, and scenarios (including fine-tuning) in the current work, which would make it a significant contribution to the community. Unfortunately, they seemed not to be present in the revision, nor in the author's immediate action plan, which prevented me from giving a higher score.

Fine-tuning

One of the concerns raised by reviewers (JL1w, SABf) was that we were not fine-tuning the scFMs in this paper. We would first like to note that even in the absence of fine-tuning, if a signal relative to perturbation effect prediction were present in the sc-embeddings, it should still be measurable in a zero-shot setting. If no performance improvement comes from foundation model embeddings and all observed performance gains are due to fine-tuning, this contradicts the guiding principle that these models learn useful biological representations during pre-training. MLP probes are commonly used as an approach to answer the question of information content of embeddings in NLP and CV, as they evaluate representation quality while removing confounding effects pertaining to task-specific prediction heads, which introduce inductive biases [1,2,3]. We also note that previous work has investigated the performance of fine-tuned versions of scGPT and scFoundation [4]. This study found that simple linear baselines still outperformed these models, suggesting that fine-tuning does not fully address their limitations. Because fine-tuning performance is addressed in other studies and because it fundamentally goes against our approach of establishing existing information content, we do not include fine-tuning in our study design. The findings we present therefore highlight that scFM embeddings have not learnt useful biological patterns pertaining to the perturbation prediction task. This in and of itself is an important finding.

[1] https://openreview.net/forum?id=SJzSgnRcKX

[2] https://arxiv.org/abs/2103.00020

[3] https://www.biorxiv.org/content/10.1101/2023.10.19.563100v1

[4] https://www.biorxiv.org/content/10.1101/2024.09.16.613342v1

In conclusion, we thank the reviewers for their thorough evaluation and constructive feedback. Through this rebuttal, we have addressed several key concerns:

1. We have clarified our rationale for evaluating models in a zero-shot setting, emphasising that this approach specifically targets the assessment of learned biological representations during pre-training. The observed limitations in perturbation effect prediction suggest that current scFMs may not be learning the intended biological principles during pre-training, an important finding for the field.
2. We have expanded our evaluation to include:
   • GEARS as an additional baseline, demonstrating that even purpose-built architectures face challenges with distribution shift Results on dual-gene perturbations, which further validate our initial findings while highlighting the increased difficulty of predicting complex perturbation effects
   • Analysis of MLP parameter counts, confirming that our findings are robust to variations in model capacity
3. We have provided a detailed justification for our dataset selection, explaining that the Norman dataset's strong perturbation signal makes it particularly suitable for evaluating whether scFMs capture perturbation-related biological information. We are actively working to extend our analysis to the Replogle dataset to further strengthen our conclusions.
4. We have defended SPECTRA's value as an evaluation framework, emphasising its ability to systematically quantify distribution shifts and provide insights that random splits cannot capture. This approach is particularly relevant for real-world applications where biological variability poses significant challenges.

Our findings suggest that instead of pursuing general-purpose scFMs, perturbation-specific design principles should be incorporated into scFM model architectures. We believe our benchmark makes a valuable contribution by providing standardised evaluation protocols that can guide such future developments. We have open-sourced these tools and the framework as a resource for the research community.

We look forward to incorporating these clarifications and additional results into our revised manuscript. We believe these changes will strengthen our paper while maintaining its core contribution: a comprehensive evaluation framework that highlights both the current limitations and future opportunities in single-cell foundation models.

SPECTRA and Dataset Splitting

We appreciate the concern regarding the impact of dataset splitting (yfBT). However, we believe SPECTRA’s contribution extends beyond what random splits or current experimental designs offer. In real-world applications like Perturb-seq, unpredictable biological distribution shifts are common. Random splits are not inherently easier or harder; a random split could align with a sparsification probability of 0.1 or 0.7, representing significantly different train-test dissimilarities. However, systematically quantifying this dissimilarity makes it possible to contextualise performance results and compare model robustness across different tasks or datasets. Random splits may therefore fail to capture these distribution shifts, often leading to overly optimistic assessments of model performance. SPECTRA addresses this by generating splits of increasing difficulty, starting with random splits and progressing to more challenging ones based on train-test similarity. This approach ensures that SPECTRA encompasses random splits but also evaluates a model’s robustness across a broader spectrum of biological variability. SPECTRA does not artificially inflate task difficulty but rather controls for the similarity between train and test split. This indeed increases test difficulty, but not artificially. While the train-test similarity control may lead to smaller datasets in harder splits, this is a reflection of the redundancy inherent to the dataset, not a limitation of SPECTRA. Moreover, SPECTRA is not just an evaluation tool but a valuable guide for experimental design. For example, in Perturb-seq, SPECTRA can reveal if a model generalises well within pathways but struggles across pathways, guiding researchers on which genes to test next. Random splits, by contrast, obscure such insights due to high train-test similarity. We argue that SPECTRA’s ability to highlight weaknesses in generalisability, simulate real-world distribution shifts, and inform experimental design makes it a high-impact addition to evaluation pipelines, especially for tasks like Perturb-seq where navigating distribution shifts is critical.

Dataset

One of the main concerns raised by reviewers (JL1w, SABf, 3E1k, yfBT) is that our results were only obtained on single-gene perturbations from the Norman dataset. We want to clarify that we originally focused on the Norman dataset because it has been shown to contain the strongest perturbation effect signal for genetic perturbations [See Figure 3.b. in [5]] and is therefore best suited to our analysis. We hypothesised that if scFM embeddings contain biologically relevant information for perturbation effect prediction, it would be more likely to be contained in a dataset with a strong perturbation signal. Positive results in this context would indicate that biologically relevant perturbation effects can be predicted accurately with zero-shot scFM embeddings, while negative results would provide greater confidence that the observed limitations reflect true negatives, rather than being due to low-quality signal inherent to the dataset tested. In line with this argument, the Replogle dataset, despite containing a greater number of perturbations (~4,000 K562; ~200 RPE1), demonstrates much diminished perturbation effects signal overall [See Figure 3.b. in [5]]. This makes it less well-suited for drawing strong conclusions about scFM embeddings in their current form. We have expanded our analysis to include all the combinatorial perturbations in the Norman dataset, where originally we only considered single-gene perturbations. We present the updated results in Tables 3 and 4.

Table 3: Perturbation effect prediction evaluation across 2000 HVGs.

Table 4: Perturbation effect prediction evaluation across the top 20 DEGs.

Note: For two-gene perturbations split 0.5, not enough perturbations passed our quality control to properly define the split.

The relative performance with respect to the mean baseline increases between single and double genes, however, the baseline gene expression MLP still performs best. This aligns with our prior observations that generalisability decreases for perturbations that are strong or atypical as indicated by higher AUSPCs. While the additional results provide further evidence of the trends observed in our initial analysis, they also underline the inherent limitations of current scFM embeddings in addressing the complexities of dual-gene perturbations and emphasise the necessity of future efforts to improve their robustness and biological relevance. Furthermore, we are currently working on including the Replogle et al. dataset in our benchmark. While memory constraints currently pose challenges to processing this dataset, we are working on solutions and aim to report these results as soon as possible.

[5] https://www.nature.com/articles/s41592-023-02144-y

MLP parameter count

Another concern raised (JL1w, 3E1k) is the discrepancy in the number of parameters of the MLPs used to probe the scFM embeddings, given that these embeddings have a variable size which is a parameter of their corresponding scFM. To alleviate this concern, we trained a range of MLPs with increasing parameter count on our raw expression data to verify the effect of parameter number on results. We report our findings in Table 2, where it can be seen the increase in parameters has no effect on the MSEs obtained. We also include additional results using scBERT and scFoundation embeddings as input, with increasing parameter count. This further strengthens our conclusion that, despite an increase in the expressiveness of the underlying MLP model, the scFM embeddings do not contain sufficient biologically relevant information from the pre-training stage to provide increased performance for this task. Despite the disparities in parameter count, our current approach maintains alignment with the standard practice of tailoring MLP architectures to the embedding size to ensure sufficient expressiveness for each model. This guarantees that no embedding is artificially constrained or unfairly advantaged due to mismatched architecture capacity.

Table 2: Train and Test set results with MLPs of increasing parameter count

Additional Baseline

To answer reviewers' concerns about including additional baselines (SABf), in Table 1 we include results where we have implemented GEARS on the raw expression, using SPECTRA splits. GEARS achieves strong predictive performance, outperforming baseline models, consistent with its purpose-built architecture specifically designed for perturbation effect prediction. Similarly to other models, GEARS achieves better performance in prediction of the overall perturbation effect compared to the effect on the top 20 differentially expressed genes (DEGs). This is in line with our previous observation that in general, fine-grained perturbation responses are harder to predict. We again observe an increase in MSE with increasing sparsification probability, a pattern consistent with our other findings. This shows that even highly specific perturbation prediction models suffer under increasing train-test dissimilarity, emphasising the difficulty of generalising to out-of-distribution perturbations. This underscores the need for robust evaluation frameworks that take distribution shift into account to assess model generalisability and highlight both current limitations and future directions for research. Alongside the finding that zero-shot scFMs fail to provide meaningful representations for perturbation prediction, these results suggest that, instead of focusing on building general-purpose scFMs, the successful prediction of perturbation effects requires incorporating perturbation-specific design principles into the architecture of scFMs.

Table 1: MSE for perturbation effect prediction with GEARS trained on raw expression data across increasing sparsification probabilities

Overall MSE

Top 20 DEGs