We thank the reviewer for the thorough review, the constructive comments and the stimulating question. We will go through these points in what follows.

Weaknesses

We thank the reviewer for the constructive comments: we propose to implement all the 4 main points suggested, as we agree they will be valuable improvements in the flow of the discussion. We add two specific remarks:

• we already provide a definition of active path in Definition 5 (L 485-486 “a path $\pi$ […] is active w.r.t. …”). For reference, the definition is taken from “Causal Reasoning with Ancestral Graphs”, 2008, Zhang; though we agree that using “a path that is not blocked” in the main text simplifies the life of the reader, as this is more commonly used.
• For the residuals estimation (eq. 12) we use kernel ridge regression, as proposed in NoGAM original paper. Note the target function in least squares regression is the expectation of the target given the covariates, i.e. $\mathbb{E}[Y|X]$ for generic X,Y variables, where Y is the target. This is consistent with our theoretical analysis, where the residuals of the score $\partial_{X_j} \log p(X)$ in equation 12 are defined as $\partial_{X_j} \log p(X) - \mathbb{E}[\partial_{X_j} \log p(X) | R_j]$ which are indeed the least squares regression residuals. We agree that the discussion about residual estimation is relevant, and we will add it in the main text of the paper.

Concerning the minor points raised, we also agree that they would make a valuable addition, we implement all points up to the last. For a title that better conveys the meaning of the paper we will consider something in the flavor of “Identifiability of latent variable additive noise models with score matching”

We are eager to mention the reviewer’s contribution in improving the clarity of our work in the acknowledgments of our paper.

Questions

This question is very interesting. We will start focusing on [1], and then will consider [3, 4], as these are the three papers we are familiar with (and where I hope I can provide a compelling answer).

1. Comment about connections to [1] (GraN-DAG): one intuition about score methods for causal discovery is that they are based on the idea that if we could access the individual kernels composing the Markov factorization, we would have all the information to recover the causal relations (assuming identifiability as in restricted ANMs). Intuitively, this is achieved as follows

   • Take the logarithm to transform the Markov factorization in a summation
   • Take the partial derivatives of this log-likelihood: partial derivatives at one time (i) isolate the kernel of leaf nodes (i.e. find the causal direction, as in our Propositions 2 and 3) (ii) inform about the dependency of each kernel from the other variables (i.e. can replace conditional independence tests, our Proposition 1).

   GraN-DAG uses the Jacobian to do the equivalent of (ii), i.e. replace conditional independence testing with sparsity of the Jacobian to find some kind of connection strength (see equation 15 of their paper). The reason why no logarithm appears is that while in our case we go through score matching estimation to try to access kernels of the markov factorization, they model such kernels directly with a neural network.

   This connection between their Equation 15 and Spantini et al., 2017 is not explicit in their work, neither in any existing paper, to the best of our knowledge, so we believe that this analysis would make a valuable addition to our paper, as a discussion in the appendix. Moreover, one interesting conjecture we can make based on these considerations is that there may be a possibility to generalize GraN-DAG to do neural causal discovery in confounded scenarios with the same guarantees that we provide. Please take this with a grain of salt.

2. Connection with [3,4]: as [3] is the nonlinear version of ICA-LiNGAM proposed in [4], we believe that connections about one work would also apply to the other. In general, the two works deal with different types of Jacobian, compared to ours or to the Jacobian found in [1]: [3, 4] deal with Jacobian of the inverse transform (linear and nonlinear), whereas our work concerns the Jacobian of the score. Despite this note, we find some point of contact may exist: Lemma A.1 of [3] proposes a condition to detect the presence of a directed path from $X_i$ to $X_j$ based on the condition $\frac{\partial f_j(X)}{\partial N_i} \neq 0.$ This is not immediately related to our Jacobian of the score, where we can identify the Markov network by checking $\frac{\partial}{\partial X_i X_j} \log p(X)$: though, a closer relation appears in the nonlinear Gaussian additive noise model case, where for a sink node we have $\frac{\partial}{\partial X_i X_j} \log p(X) = -\frac{1}{\sigma_j^2}\frac{\partial f_j(X)}{\partial X_i}$ and edges are detected by the condition $-\frac{1}{\sigma_j^2}\frac{\partial f_j(X)}{\partial X_i} \neq 0$ (note the similarity with the condition of Lemma A.1 in [3]). As they consider the partial derivative with respect to noise $N_i$, they can find direct paths, while as we consider partial derivatives with respect to $X_i$ , we can find direct edges (this statement has subtleties: the one we notice here is that we find a direct edge if and only if we know that $X_j$ is a sink, else we simply find an edge in the Markov network). Please refer to Lemma 1 in [Montagna et al., 2023a] for the derivation of the latter expression and all relevant details. We can conjecture that, in SCMs more general than nonlinear Gaussian ANM, both our method and their methods find a way to use partial derivatives of the mechanisms to probe the presence of edges in the graph. This may also be something worth adding and expanding in our paper, for a profound view on the connection between our work and existing methodologies.

We thank Reviewer TL8k for their thorough review and the valuable comments therein. One important criticism exposed by the reviewer is that the contribution of our work is limited in relation to the existing literature: in this regard, we address to general response and the first point of our response in the section below.

Weaknesses

• “The authors claim that AdaScore is able to handle a broad class of causal models (l54), but three out of four possible situations are direct applications of existing work. …“. In reply to points 1) and 2) that are made, the fact that 3 out of 4 possible situations are extension of existing works, does not mean that the 4th one is not there. Our main contribution and the relation with previous work, are largely discussed in the paper and in the rebuttal, particularly in the first point of the general response.

  Point 3) of the review says that the only novel part “may only apply to few edges.” We politely but firmly disagree with this statement: the number of edges whose direction is theoretically identifiable by AdaScore but not identified by methods such as FCI depends on the problem at hand: there is no ground to overlook this because of the subjective belief that this does not happen frequently. The reviewer agrees that there is a novel contribution w.r.t. the existing literature (from the review: “Adascore does generalize NoGAM to orient unconfounded edges of the PAG returned by FCI”;). The reviewer also acknowledge that this is not done by any other existing methods (throughout the review, it is mentioned that our method, and thus our theory, is more general than CAM-UV, FCI, and NoGAM. We mention that it is more general than RCD and LiNGAM, also). To these points we add that our work is the first one that provides identifiability theory for latent variables causal discovery with score matching, a recent and lively line of work (Montagna et al. (a, b, c), Amin et al., Zhu et al., Rolland et al.) discussed in the paper (L35) with connection to the field of causal representation learning Varici et al. (2024). The only point raised which may undermine our contribution is that our theory and method are not relevant given that the identifiability of the causal direction “may apply to few edges”, which is a purely subjective consideration, as it depends on the problem at hand and can not be simply disregarded. Given these remarks, we kindly ask the reviewer to reconsider their score on the contribution, and the global score, accordingly.

• “The authors state that they show how constraints on the Jacobian of the score can be used as conditional independence testing […]”. We do not state that the Jacobian of the score can be used as conditional independence testing or that we propose a statistical test. We say that we use the Jacobian of the score as an alternative to conditional independence testing (L9-10, L41-42): while constraint-based methods (e.g., PC) usually rely on conditional independence testing (e.g. Zhang et al., 2011), instead, we consider sparsity of the Jacobian of the score. These are two different statistical problems. Our task concerns the estimation of the Jacobian of the score (Rolland et al. for details on such estimation, Zhu et al. for its statistical efficiency, both references are in the paper), and using a t-test for finding entries of such Jacobian with zero means: the t-test is not part of our novel contributions, we will add a citation e.g. to Walpole et al., (1972) for relative details.

  Later, the reviewer writes “[…] the extent to which this is done is only by noticing the equivalence between conditional independence and the corresponding zero in the Jacobian term (noted in [1, 2])”: we don’t merely note the equivalence. We prove it as a corollary of Spantini et al., (see the 1st point of the general response). Concerning the relation with [2],referenced in our paper, they proposed a version of our Proposition 1 in the limited setting of nonlinear additive Gaussian noise models, while our result holds under generic SCM.

• “In L193, the authors state, "we remove the nonlinearity assumption (of [3]) and make the weaker hypothesis of a restricted additive noise model", but 1) this is a stronger assumption [...]” [3] assume a nonlinear restricted additive noise model, while we assume a restricted additive noise model. Allowing for linear and nonlinear mechanisms is weaker than allowing for nonlinear mechanisms only.

Limitations

We remark the adaptivity of AdaScore in the sense described in the paper, meaning that it can perform with theoretical guarantees on latent variable models, still identifying the causal direction when this is theoretically possible, which is not done by other methods (FCI, NoGAM, LiNGAM), or it is done under more restrictive assumptions (CAM-UV, RCD). We remark that this does not require any user specification. The user specification simply allows to turn off the search of latent variables, if this is thought as unnecessary. This is all described in the main paper (see L263-265, and L266 where “We […] describe the version of our algorithm whose output is a mixed graph”, i.e. AdaScore when the user is willing to account for latent variables). We could easily remove the possibility of interacting with AdaScore to provide prior knowledge, but this would be a poorer design choice (e.g. see DoDiscover for a discussion about this point from an existing library).

We renew our thanks to the reviewer for the points raised, and are willing to incorporate suggestions to better reflect the nature of our contribution, naming Proposition 1 “Corollary 1”, which better reflects its link with Spantini et al. If our response is satisfactory relative to the points raised, we kindly ask the reviewer to reconsider their score to our paper.

We thank the reviewer for the time and effort in reading our paper. One important concern unaddressed by our general response is about the limits of our experimental evaluation: we point to the first bullet in our response below, and the experimental results in the PDF of the rebuttal.

Weaknesses

• “The experimental results show limited added value according to the one metric chosen (SHD) in a synthetic setting. They are not comprehensive enough …” We highly regard this suggestions and will include real data experiments in our paper. Our results on bnlearn datasets presented in the PDF of the rebuttal (Fig. 2) show that AdaScore always outperforms CAM-UV and RCD, outperforms NoGAM on 3 out of 4 datasets, and outperforms LiNGAM on 2 out of 4 datasets. Thus AdaScore is a practical alternative with better theoretical identifiability guarantees than any of the other methods we consider. Though, please note that real and synthetic benchmarks for causal structure learning are documented as far from being satisfactory, see e.g. “Self-Compatibility: Evaluating Causal Discovery without Ground Truth”, 2024, Faller et al; experimental results should be taken with a grain of salt.

• “In line 74, the authors state that faithfulness is assumed (I believe, it is kind of hidden in the background notions).” In L74 we write “we assume that the reverse direction $X_i \perp X_j | X_Z \implies X_i \perp^d_{\mathcal{G}} X_j | X_Z$ hold”: this is the definition of the faithfulness assumption (see our references therein), so the faithfulness assumption is not hidden in background notion, but explicitly made for the first time in L74.

  ”If that is the case, the method adopts the same assumption as FCI, plus ANM. …”. We point to the second bullet of the general response for the answer on this point. Concerning the request of an experimental comparison with FCI, consider that the two methods have different identifiability guarantees and so different types of graphs as output. In this sense, FCI is not the ground for a fair comparison, as it loosens the assumptions at the price of weaker identifiability guarantees. We propose to compare the FCI F1 accuracy in skeleton recovery versus AdaScore F1 accuracy in skeleton recovery when inference is done according to Proposition 1 (at the base of the skeleton inference in AdaScore), which does not require assumptions on the SCM. Experimental results are shown in the PDF of the rebuttal (Fig. 3): we see that AdaScore consistently outperforms FCI in the task of recovery of the PAG skeleton.

• “Proposition 1 is a rather trivial application of the more general lemma in Spantini et al. and it does not specify the required faithfulness assumption to obtain the result from Eq. 6” For the first part (that Proposition 1 is a trivial application of Spantini et al.) we point to the first bullet in the general response of this rebuttal. For the second point made about Spantini et al. and the faithfulness assumptions, some clarifications are needed: Eq. 6 in the paper is the Lemma of Spantini et al. (and not obtained from it, as written in the review) and does not require faithfulness assumption. To avoid any confusion, we propose to write “Note that this result does not make use of faithfulness assumption” right after Eq. 6. Proposition 1 is a corollary of Spantini et al. that makes use of the faithfulness assumption.

• “[…] scalability of the method is not shown nor discussed.” Scalability is discussed in Appendix C.3 and Proposition 5, also referenced in the main manuscript, L 284. We propose to add experiments with larger number of variables to the camera-ready version.

• “The model used to estimate residuals is not discussed, nor the assumption that the chosen model fits the data adequately to correctly estimate residuals, and what is needed to assess this.” We use kernel ridge regression, as proposed in NoGAM original paper. Note the the target function in least squares regression is the expectation of the target given the covariates, i.e. $\mathbb{E}[Y|X]$ for generic X,Y variables, where Y is the target. This is consistent with our theoretical analysis, where the residuals of a variable $V_k$ are defined as $V_k - \mathbb{E}[V_k|V_{Z\setminus \set{k}}]$ (Equation 15) which are indeed the least squares regression residuals. We agree that this point is worth discussing, and will add it to the main text of our paper.

Questions

• “Line 28 …” As causal directions can not be identified at parity of assumptions, we do not improve upon FCI identifiability guarantees in that setting. Better identifiability comes from additional assumptions, as discussed in the general response.
• Line 44 … Yes, thank you.
• ”Line 217: ..." “Recentered” is intended in the sense that we have the random variable $f_i(V_{PA_i})$, the observed causal effects, that are shifted by $g_i(U^i)$ addition (see equation 13). We will provide a clearer formulation of this point.
• ”Line 221: ..." Although this certainly subsumes the setting of RCD, it does not necessarily. We can have a mixture of linear and non-linear dependencies between $X_i$ and its observable parents and also arbitrary relations to the unobserved parents.
• ”Line 257: ...” Our method can provide no further orientations than the ones that already follow from the Markov-equivalence class. We will write this clearer form in the paper.
• ”Line 279: ...” See section C.2
• ”Line 280 ... We also observed that pruning has a non-negligible impact. We will add the ablation study in the camera-ready version of the paper.
• ”Alg. 1: do you start from a disconnected graph? ... As we describe more explicitly in algorithm 2, we start by adding undirected edges between $X_i$ and $X_j$ if $X_i \perp X_j | X_V \setminus \{X_i, X_j\}$. We will make this more clear in the main text in the updated version.
• “Line 305: ... We will add statistical tests to the empirical evaluation in the camera ready version.

For completeness, we present the table of p-values of all the MW tests on experiments relative to dense graphs. Please notice that adding the random and fully random baselines to the experiments changed the random state of our running relative to Fig. 1 in the main paper. Results are qualitatively the same, but we observe a few mismatches due to the variance resulting from different random states.

Alternative: less

Linear fully observable model

nonlinear Fully observable model

linear Latent variables model

nonlinear Latent variables model

Alternative: greater

linear Fully observable model

nonlinear Fully observable model

linear Latent variables model

nonlinear Latent variables model

We thank the reviewer for the active participation in the discussion. We are surprised that the reviewer specifies only now that all experiments with less than ten variables and in sparse settings are not meaningful in analyzing the algorithm performance, as this point was never raised in the review or in the authors-reviewers discussion, where the requests were limited to a better analysis of scalability (without even mentioning density as a discriminating point). Regarding the claim that these are the setting common in the literature, we point to CAM-UV paper (the work closest to ours), where they have experiments on at most nine variable, and one semi synthetic experiment on ten variables. Further, we notice that the analysis of the reviewer on the goodness of AdaScore compared to other methods is not based on statistically significant statements: e.g. CAM-UV is better than AdaScore with statistical significance in one experimental setting out of 32. Yet, the reviewer concludes that CAM-UV is better than AdaScore on inference on linear latent variables models, which encompasses 8 out of 32 of our experimental settings. As the statistical tests were an explicit and well-motivated request from the reviewer to be able to make statistically significant conclusions, we will stick to the conclusions supported by the statistical tests.

In this spirit, we provide the summary table with the level of the test $0.05$, as asked by the reviewer:

These results contrast with the analysis of the superiority of certain methods made by the reviewer. Specifically, we respond to specific points made by the reviewer and not supported by the empirical evidence:

• “Linear fully observable on par with random”: we have the following list of p-values for alternative : less on linear fully observable models comparing with random: 0.0002, 6.49e-05, 0.0005, 0.31, 5.41e-6, 5.41e-6, 5.41e-6, 5.41e-6. The only case in which p-value is higher than 0.05 is 0.31. The alterantive : greater pvalues instead are never below 0.05 threshold for random, and not even close: 1, 1, 1, 1, 0.99, 0.99, 0.99, 0.71. AdaScore is not on par with random.
• “linear latent: worse than camuv”: as already discussed, there is only one time in which alternative :great in comparison with CAM-UV, in the linear latent setting, goes below 0.05. The claim is not supported by empirical evidence

Further, the reviewer writes that the linear latent case is the focus of our work: we disagree, both the linear and nonlinear settings are equally important, as the point of the paper is indeed relaxing linearity and nonlinearity assumptions on the causal mechanisms, in the context of latent variable models. If we suggested that linear latent variable cases were the focus of our work, we kindly ask the reviewer to point to the specific parts of the paper.

Finally, concerning the claim that contributions need to be clarified, we again point to L48, the paragraph starting with “our main contributions” written in bold. Concerning the claim that assumptions need to be clarified, the only specific criticism expressed on this point was relative to the positioning of faithfulness assumption, that could be more upfront. We already addressed this and remark on our offer to make a separate paragraph.