*3.2. Velocity Variations*

Linear regression analysis (see methods section) reveals a coefficient of determination R<sup>2</sup> > 0.90–0.99 for all cumulative time series (Figure 4), indicating steady surface movements with no obvious acceleration phases within the observation period. Approximately 2.6% of the analyzed time series showed an R<sup>2</sup> ranging between 0.90–0.96, however, this is primarily due to under-sampling of the time series where data points were missing (Figure 5).

**Figure 4.** Spatial distribution of coefficient of determination (R2) from the linear correlation analysis.

**Figure 5.** Spatial distribution of sampling completeness.

Extracted residuals of the time series from linear fitting revealed that the residuals of the glacier surface velocities fluctuate around the long-term (3.2 years) linear trend. Spatial correlation and cluster analysis of the residuals shows that these variations in surface velocity have a spatial and temporal variability with respect to zones in the central, eastern, and western parts of the glacier, as well as zones above 5900 m a.s.l. (Figure 6).

Surface energy and glacier mass balance models (SEMB) across the glaciated Andes [1,36,37] collectively show that the glacier response has a strong spatial variability due to its sensitivity to regional and local meteorological and topographic factors [2,37,38]. In this connection, the observed spatio-temporal variation in our glacier velocity residuals may be interpreted as being closely linked to meteorological and topographical variables, which in turn influence the variability of the glacier's surface energy, mass balance, meltwater production, and effective stress states at the glacier bed.

#### *3.3. Spatial Correlation and Cluster Analysis*

A spatial correlation and cluster analysis following the procedure described in Sections 2.3 and 2.4 were applied to elucidate details of the factors influencing the spatiotemporal velocity patterns in different parts of the glacier (Figure 7). It is important to note that this analysis is performed independently of the previous correlation evaluations and is performed directly on the fitted time series as an additional investigation using an unsupervised machine learning algorithm. From these analyses, three clusters, each with a distinct behavior, were identified from the residuals (Figure 6b). These clusters may be considered as two end-members; an intra-annual and inter-annual cluster, with a transition cluster between the two (Figure 7c), revealing a behavior similar to the results of the correlation analysis (Figure 7a,b).

**Figure 6.** (**a**) Cumulative surface velocities from six representative points, (**b**) Residuals of the linear fit displaying an intra-annual (seasonal) pattern, from three representative points above 5500 m a.s.l. (No. 143), in the East (No. 2442), and West (No. 752) of the glacier, (**c**) Residuals of the linear fit displaying an inter-annual pattern, from three representative points for the glacier response in the central, south-exposed part of the glacier (Nos. 702, 1522, 1574), (**d**) meteorological data; Temperature T, Humidity H, Precipitation P (Palcacocha weather station, altitude 4607 m.a.s.l., latitude 09◦24 09.3211"S, 77◦23 07.0258", Instituto Nacional de Investigación en Glaciares y Ecosistemas de Montaña INAIGEM), and (**e**) Sea Surface Temperature Anomalies (SSTA) shifted by 3 months (see main text) for the NOAA Niño 3.4 region.

**Figure 7.** (**a**) Map of correlation coefficients between residuals of the long-term linear fit through the velocity time series and precipitation, measured from the Palcacocha weather station, altitude 4607 m.a.s.l., latitude 09◦24 09.3211"S, 77◦23 07.0258", INAIGEM, (**b**) Map of correlation coefficients between residuals of the long-term linear fit through the velocity time series and SST, (**c**) Results of the cluster analysis show the distribution of intra, inter and transitional clusters.
