We warmly thank the reviewer for his comments.

Regarding the multi-volcano PVHA, we would like to note that with “the first multi-volcano PVHA” we referred to the novel methodology in volcanic hazard assessment which provides the aggregation of a set of volcanoes in calculating the mean annual frequency with which a certain hazardous phenomenon (in this case tephra fallout) exceeds a given threshold at a given location during a given exposure time window.

Actually, in the case of Becerril et. (2014), the authors provided a long-term multi-hazard assessment considering lava flows, pyroclastic density currents and fallout for a single volcano (even though it has a monogenetic field). In our case we considered tephra hazard only, applied to the three Neapolitan volcanoes. However, we explained better this point in the Introduction. 

We very appreciated the comments of the reviewer, and we thank him for helping us to clarify the following points:

1) We here provide an additional explanation aimed to emphasize the similarities to tsunami hazard modelling and other geohazards. In the revised manuscript we will add the following paragraph: “Other geohazards, such as for Probabilistic Tsunami Hazard Assessment (PTHA) or Probabilistic Seismic Hazard Assessment (PSHA) methodologies, provide a framework for assessing the exceedance provability of a given measure of the intensity of the phenomeno (e.g. tsunami wave height, peak ground acceleration) at a particular location within a given time window. As for volcanic eruptions, historical catalogs are usually incomplete, and thus it is usually adopted a computational hazard scheme, based on the combination of probabilistic source models and empirical or numerical models of propagation of the hazardous phenomena (e.g., Grezio et al., 2017, Gerstenberger et al 2020). Sometimes, the explicit  numerical modeling of individual scenarios is required, like in the case of tsunamis or of the volcanic eruptions(e.g., Selva et al., 2016; Grezio et al., 2017). In all cases, a common description based on the quantification of the mean annual rates of exceedance is possible, making possible an explicit comparison among the different hazard and consequent risks (Gruntal et al 2006; Marzocchi et al. 2012). Beside volcanic activity, other geohazards (as landslides, meteorological events) could be treated within the same computational framework”;

2) Regarding the rewiever's comment "While the Authors made everything possible to show the data sources clearly, the amount is so overwhelming, that the reader can easily be get lost among the various level and type of data utilized. Maybe a simple overview table could help for further clarification where the reader quickly can see the exact maps, papers, direct measurements of the utilized parameters been delivered", we would like to stress that the mean Annual Probabilities for each volcano are based on the study of Selva et al. (2022), who report a detailed catalog of eruptions over a long reference period for each volcano. Other parameters were taken from Sandri et al. (2016) for Somma-Vesuvius and Campi Flegrei, and from Primerano et al. (2021) and Selva et al. (2019) for Ischia. The above studies already presented reviews and several tables that summarize the sources needed for such pieces of information. We think adding a further composite table would not be easy and useful because it would become difficult to be presented and read;

3) Regarding the rewiever's comment "The manuscript is in its early section introducing the “eruption size class” concept. This gradually became some sort of core of the simulations, but I find not evidently clear what those classes actually mean.", we would like to stress that with the term “eruption size class” we referred to the broad range of possible eruptive sizes identified by the total erupted mass which is used to define the eruption magnitude. Following Sandri et al. (2016), we consider splitting the eruptive size range into a few classes that can be linked to representative members like the classical approach used in past studies. These classes ideally span the general range continuously, whereas representative members, by definition, discretize it. We added a short explanation in the revised manuscript;

4) Regarding the rewiever's comment "I was wondering if a narrative table would not be a good addition where the reader (basically those approaching to the manuscript or try to understand it from more traditional geology background) could see clearly what those values meant.", as explained above, we think that such a kind of very heterogeneous table would not be feasible and useful, and in any case, all the relevant references are clearly indicated;

5) Regarding the rewiever's comment "the manuscript has also few very specific terms unlikely would make sense for those readers not deeply involved in geomathematics. Terms like “intra-size-class aleatory variability” would be worth to explain for clarity.", as explained in the previous points, we split the eruptive size range into a few representative classes. Within each class there is a variability in size that is indicated as “intra-size-class aleatory variability”, as described in Sandri et al. (2016). In this regard, the “intra-size-class” variability represents the aleatory variability due to combinations of parameters characterizing eruptions which belong to the same eruptive class.  As for the previous point, in the revised version, we will add an explanation about this point;

6) Regarding the rewiever's comment "One of the main selling points of the manuscript is to combine the three main hazard source fallout hazards and develop a united hazard model. Somehow, I see why this is important, but I think a clear paragraph to show the benefit of doing this would really set the manuscript far better. ", we agree with the reviewer, therefore we will provide the following paragraph in the revised manuscript:  “The total mean hazard maps shown in Fig. 3d provide a unified representation of the potential tephra ground load coming from one of the three volcanoes. This means that, for different ARPs, we can evaluate the main tephra ground load exceeding the defined thresholds, in each point of the whole domain. As expected, the shorter ARP (100 years) represents the most frequent hazardous scenario related with the most frequent eruptions having lower magnitudes, while for the longer ARP (1000 years) prevails heavier tephra load distribution related with less frequent eruptions having the highest magnitude. In this way, such overall tephra hazard picture is more useful for mitigation actions (e.g., road viability, roofs damaging) with respect to the previous tephra hazard maps which accounted for a single volcano.”;

7) Regarding the rewiever's comment "According to the info on the geological map, Ischia also had small or medium tephra providing eruptions in the same timescale than the other two volcanoes, so just not really clear why then they are not included. Really, just an explanation line would clear this.", we thank the reviewer for this comments that lets us to clarify the choice made for Ischia. We only considered the Large eruptive size class (i.e., Cretaio eruption) because it is the best characterized eruption in terms of ESP (Primerano et al., 2021) and is the only eruptive scenario which can have a significative impact beyond the island (Selva et al., 2019). On the contrary, the other eruptions having lower magnitude do not have a significant impact beyond the island (Selva et al., 2019). We will add this explanation in the revised manuscript;

8) Regarding the rewiever's comment "The hazard maps are very informative. I was wondering a lot about the concept of hazard disaggregation and the claimed similarity with the seismic source differentiation. I think a bit of explanation is needed here.", we would like to stress that, as reported in Bazzurro and Cornell (1999), the hazard disaggregation scheme permits to postprocess hazard results to display the relative contributions to the hazard of the different source. This is not trivial, as the contribution of each individual source depend on many factors, like the size and the position of the event, the annal rates of each specific size and position, and the propagation. These are all ingredients of the hazard, but their combination and their balance is not trivial at all. For example, is more impacting in one unlikely big event or a more likely smaller event with favorable wind? Is more impacting one unlikely event in a location usually upwind, or a more likely event in a location usually downwind, but with a rare variable wind direction? To answer these questions, the only solution is to post-process the hazard combination and balance the different contribution. We will try  to add a few sentences in the text to better explain this very important point;

9) Regarding the rewiever's comment "The hazard maps look good, but I would add a line to the caption with a plain language explanation. ", we would like to note that hazard maps showed in the Figures are already a graphical representation of the tephra hazard. We think that the captions are already written in a clear way and trying to make the explanation plainer may introduce misunderstandings;

Finally, we thank the reviewer for his appreciation of our work.