Determining the eruption frequency-Magnitude (f-M) relationship for data-limited volcanoes is challenging since it requires a comprehensive eruption record of the past eruptive activity. This is the case for Melimoyu, a long-dormant and data-limited volcano in the Southern Volcanic Zone (SVZ) in Chile with only two confirmed Holocene eruptions (VEI 5). To supplement the eruption records, we identified analogue volcanoes for Melimoyu (i.e., volcanoes that behave similarly and are identified through shared characteristics) using a quantitative and objective approach. Firstly, we compiled a global database containing 181 variables describing the eruptive history, tectonic setting, rock composition, and morphology of 1,428 volcanoes. This database was filtered primarily based on data availability into an input dataset comprising 37 numerical variables for 438 subduction zone volcanoes. Then, we applied Agglomerative Nesting, a bottom-up hierarchical clustering algorithm on three datasets derived from the input dataset: 1) raw data, 2) output from a Principal Component Analysis, and 3) weighted data tuned to minimise the dispersion in the absolute probability per VEI. Lastly, we identified the best set of analogues by analysing the dispersion in the absolute probability per VEI and applying a set of criteria deemed important by the local geological service, SERNAGEOMIN, and VB. Our analysis shows that the raw data generate a low dispersion and the highest number of analogues (n = 20). More than half of these analogues are in the SVZ, suggesting that the tectonic setting plays a key role in the clustering analysis. The eruption f-M relationship modelled from the analogue's eruption data shows that if Melimoyu has an eruption, there is a 49% probability (50th percentile) of it being VEI=4. Meanwhile, the annual absolute probability of a VEI=1, VEI 2, VEI 3, VEI 4, and VEI=5 eruption at Melimoyu is 4.82 x 10(-4), 1.2 x 10(-3), 1.45 x 10(-4), 9.77 x 10(-4), and 8.3 x 10(-4) (50th percentile), respectively. Our work shows the importance of using numerical variables to capture the variability across volcanoes and combining quantitative approaches with expert knowledge to assess the suitability of potential analogues. Additionally, this approach allows identifying groups of analogues and can be easily applied to other cases using numerical variables from the global database. Future work will use the analogues to populate an event tree and define eruption source parameters for modelling volcanic hazards at Melimoyu.