Aftershock Regions of Aleutian‐Alaska Megathrust Earthquakes, 1938–2021

With five earthquakes (1938, 1946, 1957, 1964, and 1965), the Aleutian‐Alaska subduction plate boundary ruptured a length of 3,548 km. We revisit these five earthquakes—first studied in detail by Sykes (1971, https://doi.org/10.1029/JB076i032p08021)—through probabilistic relocation of carefully selected mainshocks and aftershocks. Our final catalog of 324 events is established from a set of 12 mainshocks that includes all Mw ≥ 7.7 megathrust earthquakes. Using the relocated catalog, we create revised aftershock regions delimited both parallel and normal to the trench. These aftershock regions exhibit significant differences from previous studies, with the following basic findings: the 10 November 1938 Mw 8.3 earthquake extended further west, to the Shumagin Islands, and further east, into the Kodiak region, relative to the prevailing aftershock region established by McCann et al. (1979, https://doi.org/10.1007/978-3-0348-6430-5_2). The 01 April 1946 Mw 8.6 sequence was anomalously concentrated near the trench, which implies near‐trench coseismic slip that contributed to the exceptionally large tsunami. The 09 March 1957 Mw 8.6 aftershocks spanned a 1,230 km length with numerous aftershocks within the outer‐rise region of the incoming Pacific plate. The 28 March 1964 Mw 9.2 aftershocks extended east into the Pamplona thrust system (south of Icy Bay, Alaska), suggesting coseismic rupture into this region; this is consistent with coseismic static displacements, as well as current estimates of interseismic locking. The post‐1965 events we examine are all smaller than the earlier events, but since they occurred in an era of improved data collection (seismic, geodetic, and tsunami), they provide a better opportunity for assessing the link between the distribution of aftershocks and the occurrence of coseismic, postseismic, and interseismic slip.

Clearly there is demonstrated value in representing the aftershock regions of historical large earthquakes, as they provide insights into the coseismic rupture. For example, for more recent earthquakes in the same regions, one can compare their aftershock regions with coseismic slip models derived from globally and regionally recorded seismic and tsunami waveforms. The aftershock regions can also be compared with plate coupling estimates derived from modern GPS measurements and with subsurface structural heterogeneity, either on the upper plate or downgoing plate.
The underlying information used to determine the aftershock regions are hypocenters estimated from arrival time data. The hypocenters guiding the widely published aftershock regions were from , either estimated in that study or obtained from catalogs at that time. The primary purpose of our study is to apply modern methods of relocation and uncertainty estimation to aftershocks of all M w ≥ 7.7 earthquakes occurring on the Aleutian-Alaska megathrust. This results in 12 earthquakes between 1938 and 2021, six of which (1938, 1946, 1948, 1957, 1964, and 1965) were analyzed in , with the remaining events (1979, 1986, 1996, 2003, 2020, and 2021) occurring after its publication. Sykes (1971, p. 8026) explained his event selection as follows: "Events were taken to be aftershocks if they occurred within a few months and within a few hundred kilometers of the mainshock." We devote considerable effort toward documenting our selection of events, starting with 553 earthquakes that occurred within ±100 days of the 12 mainshocks. In establishing our final catalog of 324 events, we exclude events that most likely did not occur on the subduction interface (Section 4). Among the 324 events, only 70 were relocated in .
A secondary purpose of our study is to provide updated aftershock regions based on the megathrust aftershock catalog of 324 events. The previously published aftershock regions evolved in appearance from  to McCann et al. (1979) to   (Figure 1), but no details were provided as to how the regions were drawn. Sykes (1971, p. 8030) stated, "Aftershock areas are denoted by the hatched areas in Figure 4." In our study we determine new aftershock regions based on the new catalog of 324 events.
We find that the six pre-1970 aftershock regions are significantly different from those used in the literature. The most notable changes are associated with the 1964 and 1938 earthquakes. The revised 1964 aftershock region extends more than 100 km further east, into an area that experienced coseismic uplift (Plafker, 1965) and exhibits strong coupling between the downgoing Pacific/Yakutat plate and the overriding North America plate (Elliott et al., 2013). The revised aftershock region is similar to the estimated region of interseismic locking (Suito & Freymueller, 2009). The revised 1938 aftershock region extends more than 100 km east and 100 km west of the region displayed in McCann et al. (1979). Previous studies that have used the 1938 region-or the more recent 2020 and 2021 coseismic rupture regions-for understanding active tectonic structures may wish to reconsider their interpretations in the context of the revised 1938 region. Additional differences between the old and new regions are discussed and summarized at the end of the paper.  . (b) McCann et al. (1979). (c) . (d) This study, displaying the largest 5 events out of 12. Note that (a-c) also include earthquakes along the Queen Charlotte plate boundary. Figure numbers for each aftershock region in our study are listed in Table 5.
We also use source mechanisms published in papers by William Stauder and colleagues. Stauder played a leading role in the 1960s in determining focal mechanisms for globally recorded earthquakes. His primary objective was to use S waves, in addition to P waves, to improve estimates of source mechanisms. He published global compilations in Stauder and Bollinger (1964); Stauder and Bollinger (1966b) and Aleutian compilations focusing on the 1957 aftershocks (Stauder & Udias, 1963), the 1964 aftershocks (Stauder & Bollinger, 1966a), the 1965 aftershocks (Stauder, 1968a), and the tensional mechanisms of events near the Aleutian trench (Stauder, 1968b). Figure S1 in Supporting Information S1 shows all 73 mechanisms published within the four Aleutian compilations. Given the data quality and available methods of the era, we assume that uncertainties in the source mechanisms for the Stauder catalog are much larger than those for the GCMT catalog.

Mainshock Event Selection
Our selection of mainshock events is as follows: 1. We consider all earthquakes in our seismic catalog starting in 1938. Since the start of the instrumental era in about 1897, there are no definitive M w ≥ 8 megathrust earthquakes up until 1938 (see Section 6.6). 2. We sketch a boundary region that encompasses the extent of the Aleutian-Alaska megathrust as defined by Slab 2.0 (Hayes et al., 2018). The boundary region contains a buffer zone to account for uncertainties of epicenters, which is especially needed for historical events. 3. We select all events M w ≥ 7.7 within the boundary region between 01 January 1938 and 01 January 2022. This results in 18 candidate events. 4. We exclude certain events. An earthquake near Umnak Island on 02 July 1965 is assigned a magnitude of M w 7.8 in the ISC-GEM catalog, based on a seismic moment estimate in Wyss (1970). We reassign M s 6.5 from Abe (1981)-based on surface wave amplitudes-and note that the earthquake does not appear in the M > 7 compilation of Pacheco and Sykes (1992). We exclude this earthquake on the basis of its smaller size. We exclude 7 earthquakes that did not occur on the subduction interface:  Table 1).
The 1948 earthquake occurred west of the Shumagin Islands. Davies et al. (1981, p. 3825-6) wrote: "  relocated only two aftershocks of the 1948 event; the computed depths from pP-P are well constrained (again, if pP-P times reported by the ISC can be taken at face value) at 54 and 86 km. Hence, the depths and locations indicate that the 1948 shock probably involved rupture at depths below about 40 km. Thus, that event may well not have been an interplate earthquake but rather a deeper shock within the downgoing plate." Estabrook et al. (1994) reexamined the 1948 event and estimated a source mechanism based on P waveforms. They concluded that the earthquake "occurred on a shallow dipping thrust fault with a depth of about 31 km, not 60 km as originally suggested." Abe (1981) listed M s 7.5 for the 1948 earthquake. Estabrook et al. (1994) estimated M w 7.1 based on P and S waveform modeling, and they acknowledged that "this difference in magnitude may be caused by the band-limited nature of the P and SH wave recordings in 1948, since nearly all of the instruments have natural periods of about 10 s, while the duration for the 1948 event is nearly 26 s." The ISC-GEM catalog uses the 7.1 magnitude, but we adhere to Abe's surface-wave-based magnitude (7.5) in our lists and figures.
The 28 February 1979 M w 7.5 earthquake occurred in the tectonically complex region containing the Wrangell-St. Elias Mountains. Several studies examined this earthquake in its aftermath Lahr et al., 1980;McCann et al., 1980;Perez & Jacob, 1980;Stover et al., 1980  the event and concluded: "The locations and source mechanisms of these subevents and locations of aftershocks define a shallow dipping surface at the eastern edge of the Pacific plate." Their seismic moment estimates from body waves and surface waves were 9.4 × 10 19 N-m and 13 × 10 19 N-m, respectively, corresponding to M w 7.2 (7.248) and M w 7.3 (7.342).

Relocation of 553 Earthquakes Using NonLinLoc
We consider all events in our seismic catalog that occurred within our boundary region and within ±100 days of the 12 mainshocks in Table 1. At the outset, we included events prior to the mainshock because previous studies suggest that slip on the subduction interface may initiate in the days, months, or years prior to a major earthquake Ruiz et al., 2014;Schurr et al., 2014;Yokota & Koketsu, 2015). The resulting number of earthquakes is 553.
We use the robust, probabilistic, grid-search location algorithm NonLinLoc (NLL) (Lomax et al., 2000(Lomax et al., , 2014 to relocate hypocenters for all 553 events. Arrival times of most events were available from the ISC database. For 50 events, we obtained arrival times from personnel at ISC. NLL uses an efficient global sampling algorithm to map the posterior probability density function (PDF) (Tarantola, 2005) in 3D space for the hypocenter location.
The location PDF provides a complete description of likely hypocenter locations and includes comprehensive uncertainty information. A PDF solution is particularly important for earlier instrumental locations and remote areas, for which there may be few stations, poorly distributed with respect to the epicenter, and large observational errors, leading to large and sometimes highly irregular epicenter and depth uncertainties.
For each relocated earthquake, NLL thoroughly explores and maps over a geographic region the PDF for the hypocenter parameters, longitude, latitude, and depth, using tens of thousands of evaluations of a misfit function. This exploration is guided by an efficient, cascading grid-search, importance-sampling method called oct-tree (Lomax et al., 2014). To aid in statistical calculations and for ease of visualization, a compact set of about 1000 posterior samples is drawn stochastically from the oct-tree PDF grid for each event. These samples form a 3D cloud of points which can be used to represent our confidence in the solution; for example, the volume containing the densest 90% of the posterior samples is the 90% confidence region for the location. The NLL misfit function we use is based on the equal differential time (EDT) formulation (Font et al., 2004;Zhou, 1994), which is more robust than standard location methods in the presence of outliers in the data. (An outlier observation has a residual much larger than its nominal uncertainty.) Travel times for first arrival P (P, Pg, p, Pn, Pdiff, PKP, and PKIKP) and S (S, Sg, s, Sn, Sdiff, SKS, and SKIKS) phases are calculated using the TauP Toolkit (Crotwell et al., 1999) with the spherical Earth model ak135 (Kennett et al., 1995); uncertainties for all travel times are set to 2 s. The ISC catalog does not list pick uncertainties, so considering possible timing error on older, analog seismograms and the P arrival residuals remaining after modern, teleseismic relocation using station corrections and 3D earth models (e.g., Simmons et al., 2012) we set 1-sigma observation uncertainties to 1 s for P readings and to 4 s for S readings, considering their longer period and immersion in the P wave coda. We primarily rely on the insensitivity of the NLL EDT formulation to outliers to effectively ignore observations with large error. We performed the NLL EDT search over the volume bounded by longitude 160°-235°, latitude 45°-70°, and from the surface to 200 km depth. To help improve relative event locations throughout the study area, we develop static corrections for P and S readings at each station. We perform an initial relocation with NLL EDT with uncorrected readings from which corrections are obtained for each station from the mean of P and of S residuals of ≤3.0 s, for relocated events with a 68% error-ellipsoid principle axis half-width ≤50 km, root mean square of residuals ≤3.0 s, number of P and S readings ≥30, and azimuth gap ≤180°. A second pass of relocations using these static corrections form the final relocated catalog presented and discussed in this study. For more details see application of NLL using historical data (Lomax, 2005(Lomax, , 2008Silwal et al., 2018;Tape et al., 2017Tape et al., , 2021. A comparison of the differences between the NLL hypocenters and the original catalog hypocenters is shown in Figure S2 in Supporting Information S1. The events with the largest discrepancies are labeled; these events tend to also have the highest uncertainty estimates. The general epicenter differences of 10-30 km agree with estimates of absolute teleseismic location errors (Wyss et al., 2011), and some bias in epicenter is to be expected due to differences in station corrections, location algorithm, and effective weighing of heterogeneously distributed 7 of 28 stations (e.g., generally higher density of stations to the north in Europe for earlier events (Tape et al., 2017), and of nearby stations to the northeast in Alaska for later events).

Selection of Events Used to Define the Aftershock Regions
The selection criteria for events were based on NLL hypocenters as well as the following information: 1. The earthquake source mechanisms-for events in the GCMT catalog (post-1976) or the papers by Stauder and others, 2. The subduction interface Slab 2.0 (Hayes et al., 2018), and 3. The mainshock magnitudes (Table 1), which provide a proxy for slip area via magnitude-area empirical scaling (Allen & Hayes, 2017).
As described in Section S1 in Supporting Information S1, we allowed for reselection of the maximum likelihood hypocenter based on the relative location of the subduction interface; this affected only seven events, two of which were included in the final set of 324 events ( Figures S3 and S4 in Supporting Information S1).

Selection Steps
The full set of 553 events is shown in Figure 3. Next, we explain the steps to select the final set of 324 events that we assume to have occurred on the subduction interface in the aftermath of the 12 mainshocks. Table 2 lists the number of events removed at each step.
1. Remove events whose source mechanisms are farthest from being shallow thrust. The GCMT catalog reveals a variety of source mechanisms along the subduction plate boundary ( Figure 2). Among the 553 events, 273 have source mechanisms from either Stauder (46 events: Figure S1 in Supporting Information S1) or from GCMT (227 events: Figure 2). For these events, we exclude any source mechanisms having a tension-axis with plunge angle less than 25° or with azimuth toward the southeast (50° − 230°) ( Figure S5 in Supporting Information S1). This step removed 66 events, which are shown in Figure 4. 2. Remove foreshock events. Lacking source mechanisms for most of the pre-1976 events, we cannot be sure that the foreshocks occurred on the subduction interface, and therefore we decided to remove all foreshocks Figure 3. All 553 events relocated in this study and colored by depth. These events occurred within ±100 days of the 12 mainshocks (Table 1) and within the region denoted by the dashed line. Beach ball symbols denote source mechanisms from the GCMT catalog (post-1976) or from Stauder (pre-1976); solid-colored circles denote pre-1976 events that were not analyzed by Stauder. 8 of 28 from the set of 553 events. This step removed an additional 71 events, which are shown in Figure 5. We discuss foreshocks in Section 6.5. 3. Remove events that are far from the mainshock. Next, we considered how far each possible aftershock was from the mainshock. The event was removed if it was further than some number of magnitude-scaled fault lengths from the mainshock epicenter. The nominal fault length for each earthquake was estimated using the empirical relationship of Allen and Hayes (2017): l(M) = 10 −2.90+0.63M , where l is in km and the mainshock magnitude M is from Table 1. For each mainshock, we assigned a selection radius of l, 2l, or 3l. This selection criterion is shown in Figure S6 in Supporting Information S1 for the 07 May 1986 M w 7.9 earthquake. This step removed 62 events, which are shown in Figure 6. 4. Remove events far from the subduction interface. From the remaining events, we evaluated the subduction interface depth at each epicenter. If the epicenter was outside the spatial extent of the subduction interface or if the hypocenter was greater than 40 km vertical distance from the subduction interface, we excluded the event.
We also excluded events deeper than 100 km. Figure S7 in Supporting Information S1 depicts this selection criterion. This step removed 21 events, which are shown in Figure 7. 5. Remove events manually. This step removed 9 events, which are shown in Figure S8 in Supporting Information S1. Three events are 1964 aftershocks that occurred in the Cook Inlet region, about 50 km northwest of the other aftershocks. Based on the prevailing mechanisms and depths of events in this region ( Figure 2), we assume that these 3 events occurred in the downgoing slab. We excluded 4 events that occurred in the time between the non-megathrust earthquake on 19 October 2020 M w 7.6 earthquake and 100 days after 22 July 2020 M w 7.8. We excluded a 1957 aftershock (ISC 886310) based on its location near the trench, where other aftershocks slightly further south were excluded ( Figure 7). Also excluded was a M w 6.9 intraslab event on 21 October 2021. 6. Add events manually. There were 5 events ( Figure S9 in Supporting Information S1) that would have been discarded in the steps above, but we instead chose to keep them in the analysis. All of these events fall outside the extent of the subduction interface defined by Slab 2.0 (Hayes et al., 2018). One of these events (20 April 1979 M w 5.4) is the only aftershock of 28 February 1979 M w 7.5 in our analysis. The other four occur in a structurally complex region containing a set of thrust faults within the Yakutat terrane. The largest of the four events occurred on 17 May 1964 with an ISC-GEM magnitude of M w 6.1. (The published value of 5.1 in Stauder and Bollinger (1966a) is possibly a typo). The source mechanism estimated by Stauder and Bollinger (1966a) exhibits low-angle thrusting directed slightly east of north. By including these events, we suggest that it is possible that the 1964 earthquake mainshock ruptured into the complex region of thrust faults (Section 6.4).
The event selection steps are summarized in Table 2 and Figures 4-7. The full set of 324 events kept is shown in Figure 8; the full set of 229 events cut is shown in Figure S10 in Supporting Information S1. Figures S11-S13 in Supporting Information S1 summarize the distribution of aftershocks in longitude, magnitude, depth, and time. Table 3 shows the breakdown of the 324 events by mainshock. It also shows the number of these events that appeared in the Stauder compilation and within the list of 249 events in Table 1 of . Among the 46 events in the Stauder compilation that match one of the 553 events, 3 are within the final list of 324 (Table 3).

Removed Cumulative
Starting set of events -553 Step 1: remove non-thrust focal mechanisms 66 487 Step 2: remove foreshocks 71 416 Step 3: remove events far from mainshock 62 354 Step 4: remove events far from subduction interface 21 333 Step 5: manual removal 9 324 Note. Note that if the order of the steps were changed, then most of the entries would change but not the final number of 324 kept events. The relocations performed by  focused on the 09 March 1957 M w 8.6 earthquake ( Figure S14 in Supporting Information S1): 1957 aftershocks occupy 55% (137/249) of his Table 1, which spans from 04 May 1923 to 24 June 1970 (Tape, 2022). Most of his 1957 aftershocks do not appear in the ISC-GEM catalog and presumably do not meet the catalog's criteria for relocation quality.
The list of Sykes is missing 19 events (53 vs. 34) that we examined and which occurred after 01 April 1957 (Table 3). For aftershocks of the 1964 and 1965 earthquakes,  used locations from the U.S. Coast . Events (n = 66) cut because their source mechanisms were not north-dipping thrusts. Solid lines depict 5 aftershock regions from this study (1964, 1938, 1946, 1957, and 1965). Dashed lines depict aftershock regions from .

Figure 5.
Possible interplate megathrust events (n = 71) that were cut due to occurrences within 100 days prior to one of the 12 mainshocks. Each source mechanism has a thrust component, because other mechanisms were removed in the previous step ( Figure 4). and Geodetic Survey; therefore they do not appear in his Table 1. For the six mainshocks common to our study and , we relocated 244 earthquakes, 70 of which were relocated by Sykes.

Aftershock Regions for the 12 Earthquakes
The relocated catalog of 324 events (including 12 mainshocks) and their probability clouds provide the most complete picture of the spatial coverage of aftershocks of major earthquakes along the Aleutian-Alaska megathrust. This information is provided as a supplemental data set in Lomax (2022).  As discussed in the Introduction, it has proven useful to depict a set of aftershocks as a closed curve defining an aftershock region (McCann et al., 1979;). Using our revised catalog of aftershocks, we determined new aftershock regions by starting with polygon outlines containing the maximum likelihood epicenters ( Figure  S15 in Supporting Information S1). Using the polygons as a guide, we manually discretize aftershock regions in order to lessen the sharp edges of the polygon and to encompass the PDF clouds associated with the maximum likelihood epicenters. We provide the aftershock regions as digital files of longitude-latitude points (Tape & Lomax, 2022). For 6 of the 12 earthquakes (ordered from east to west), we had choices to make in creating the aftershock regions, as documented next: 1. 1979 ( Figure S16 in Supporting Information S1). This M w 7.5 event has only one aftershock in the GCMT catalog. The relocated aftershock occurred about 50 km southeast of the mainshock epicenter; Estabrook et al. (1992, p. 6589) wrote: "The southern part of the aftershock zone, while containing many aftershocks, appears not to have ruptured coseismically, but may have failed later by aseismic creep as seen in geodetic measurements." Based on the limited number of events, we did not draw an aftershock region for 1979. 2. 1964 ( Figure 10). Figure 9a shows seismic activity over the past 45 years, while Figure 9b shows the activity-including some source mechanisms (Stauder & Bollinger, 1966a)-within ±100 days of the mainshock. In drawing the aftershock region in Figure 10, we ignored one particular event-at the southwestern tip of the Kenai Peninsula-that was more than 100 km away from all other aftershocks. This M w 5.4 event (ISC 869868) occurred 6 hr after the mainshock and has a maximum likelihood depth of 29 km, which is consistent with the location of the subduction interface. Unfortunately, the event was not analyzed by Stauder and Bollinger (1966a), so we do not know its source mechanism. The only M w ≥ 5.45 event from the GCMT catalog in this region is an intraslab event (05 May 1999 M w 5.7; Figure 9a). In drawing the aftershock region, we included a somewhat isolated event under southwestern Kodiak. This M w 5.7 earthquake (ISC 868717) occurred 16 days after the mainshock. Including this event pushes the Note. Also shown are the number of our events listed in Stauder ( Figure  S1 in Supporting Information S1) and Sykes (1971, Table 1). The study of  likely used many of the same 1964 and 1965 aftershocks considered in our study; however, he did not relocate those aftershocks and therefore they were not listed in his Table 1.

Table 3
Partitioning of the 324 Kept Events, by Mainshock boundary of the aftershock region northerly, corresponding to an area overlying the subduction interface that is about 10 km deeper. This region that has produced only intraslab events in the GCMT catalog (Figure 9a). North of the mainshock epicenter, two normal-fault mechanisms from Stauder and Bollinger (1966a) result in the events being excluded. This suggests that the other 1964 events nearby, which do not have source mechanisms, could be non-interface events. We also note the only two M w ≥ 5.45 events in the GCMT catalog from this region have normal-fault mechanisms (Figure 9a). 3. 1938 ( Figure 11). The 1938 aftershock region was drawn based on nine aftershocks and the mainshock. The westernmost aftershock (ISC 902795) occurred under Unga Island 5 days after the mainshock. The event has M w 5.9 and is in the ISC-GEM catalog. The easternmost aftershocks (ISC 902776 and 902860) occurred under southwest Kodiak Island; one occurred 2 hr after the mainshock, the other occurred 28 days after. 4. 2021 ( Figure S17 in Supporting Information S1). 5. 2020 ( Figure S18 in Supporting Information S1). 6. 1948 ( Figure S19 in Supporting Information S1). 7. 1946 ( Figure 12). The westernmost event is isolated from the rest and was ignored in drawing the aftershock region. The event (ISC 898342) occurred 6 days after the mainshock and has an unknown magnitude. Our estimated epicenter (longitude −166.11) is west of the one determined by López and Okal (2006) (longitude −165.38) and west of the one in the ISC-GEM supplemental catalog (longitude −164.77). The easternmost event (ISC 898316), also listed in the ISC-GEM supplemental catalog, occurred 1 hr after the mainshock. 8. 1957 (Figure 14). Both the GCMT catalog and the events analyzed by Stauder demonstrate the preponderance of normal-fault source mechanisms near the trench and south of the trench in the region of the 1957 earthquake ( Figure 13). These events relocated just south of the trench and were excluded by our selection criteria (ISC 886310 was manually excluded). Many of the source mechanisms estimated by Stauder and Udias (1963) exhibit a considerable strike-slip component ( Figure 14). (Note that the 1963 paper does not display the source mechanisms as beach balls, and we are not aware of these being plotted in any subsequent publication.) We suspect that these aftershocks occurred on the subduction interface-therefore as shallow thrust events-and that the inferred mechanisms have some systematic errors, either on account of data errors or modeling errors. 9. 1986 ( Figure S20 in Supporting Information S1). The westernmost event is somewhat isolated from the rest. This M w 5.3 event (GCMT B051486C) has a shallow thrust mechanism and occurred 6 days after the mainshock. 10. 1996 ( Figure S21 in Supporting Information S1). 11. 1965 (Figure 15). 12. 2003 ( Figure S22 in Supporting Information S1).
We calculate the area of each aftershock region. We also calculate along-trench distances by first projecting all aftershock region boundary points onto the trench, using the digital file in Bird (2003). The along-trench length of the aftershock region is then given by the difference between the easternmost projected point and the westernmost projected point.
The differences between our aftershock regions and those previously published are best viewed in Figure 1. These differences can also be quantified, as shown in Tables 4 and S1 in Supporting Information S1. Table 4 lists areas and lengths of five earthquakes, derived from aftershock regions in three different studies. It also lists the Jaccard index between the aftershock regions in this study and those in . The Jaccard index measures the amount of overlapping area as the intersecting area of the two shapes, divided by the union area of the shapes (J = P 1 ∩ P 2 / P 1 ∪ P 2 ; see Figure S23 in Supporting Information S1). Two perfectly overlapping shapes would have J = 1; two nonoverlapping shapes would have J = 0. All along-trench distances are calculated based on the digital file in Bird (2003).
For the study of , we have calculated lengths and areas based on the solid-filled regions, not on the expanded dashed-line regions that were used for 1946 and 1938 (Figure 1a). For the study of , we have calculated lengths and areas based on the solid-filled regions, while recognizing that two of the regions (1946 and 1957) included dashed-line portions (Figure 1c).
For a pair of aftershock regions for the same earthquake, Table S1 in Supporting Information S1 lists: (a) the along-trench differences in position for the western and eastern boundaries, (b) the along-trench differences in length, (c) the ratio of areas, and (d) the Jaccard index. Table S1 in Supporting Information S1 provides details that may not be immediately apparent from Figure 1 (Table S1 in Supporting Information S1). The differences between the positions of the endpoints vary from 6 to 145 km, and the area ratios vary from 0.52 (1965) to 1.33 (1946). Table S1 in Supporting Information S1 we also compare our regions with those of , and we compare  with . The latter comparison reveals notable differences that are not addressed in the literature, notably within McCann et al. (1979), which is the primary source of the boundaries in . For example, the 1965 aftershock region used in McCann et al. (1979) and  is 138 km longer and 1.5 times the area than the one in . Presenting the aftershock regions, McCann et al. (1979McCann et al. ( , p. 1100 wrote: " Figure 5 shows the rupture areas of great earthquakes along the Alaska-Aleutian seismic zone since about 1920." No discussion was included as to why their aftershock regions differed from those of Sykes (1971).

An Aftershock Region Is Not the Same as the Coseismic Rupture Area
We use the term "aftershock region" for something that has previously been referred to as aftershock zone , aftershock area , rupture zone , and rupture area (McCann et al., 1979;. The choice of the second word-region, zone, area-is not nearly as significant as the choice of the first word-aftershock or rupture.  used the terms aftershock zone and rupture zone interchangeably, but he also identified the distinction: "Precise locations of aftershocks provide a valuable indirect method for estimating the extent of primary faulting in earthquakes (Benioff, 1955;Wilson, 1936)" (p. 8022). Wilson (1936) examined the 20 December 1932 M s 7.2 Cedar Mountain (Nevada) earthquake and concluded: "It seems logical that the distribution of the smaller shocks accompanying a major earthquake might well indicate the extent of the focal region of the main shock" (p. 1087 Note. The Jaccard index measures the amount of overlap between two aftershock regions, and it is calculated for TL and S81 (see also Table S1 in Supporting Information S1).

Table 4 Areas and Trench-Parallel Lengths of Aftershock Regions From Three Studies-Tape and Lomax (TL), Sykes et al. (1981) (S81), and Sykes (1971) (S71)-And for Five Earthquakes
catalog of Richter (1955) for the 21 July 1952 M w 7.3 Kern County (California) earthquake, Benioff (1955) wrote: "Assuming that the fault segment which was active in the production of the principal shock is effectively defined by the distribution of aftershock epicenters…" (p. 199). Thus, it has been accepted for many decades that aftershock extent provides a general picture of the mainshock rupture area. However, it was not until the 1990s that seismic data quality and coverage enabled the estimation of details of the large earthquake ruptures, such as where on the fault the greatest slip occurred (e.g., Hartzell & Heaton, 1983;Wald & Heaton, 1994).
Even with modern seismic instrumentation and modern methods, it is difficult to estimate the time-and space-dependent coseismic slip for a large earthquake, especially using only teleseismic seismic P and S waveforms (Ide et al., 2007;Mai et al., 2016;Shao & Ji, 2012). The use of additional data-such as surface deformation (from GPS or InSAR) or tsunami waveforms-helps reduce uncertainties for coseismic slip models for large subduction zone earthquakes (e.g., Sladen et al., 2010;Yue et al., 2014). While there are uncertainties in the locations of high-slip patches of large earthquakes, there are also uncertainties in the locations of aftershocks. Noting these uncertainties, Wetzler et al. (2017); Wetzler et al. (2018) examined 101 M w ≥ 7.0 subduction zone earthquakes between 1990 and 2016 to explore the relationship between aftershock locations and areas of high coseismic slip. They identified a low occurrence of aftershocks in the regions of highest coseismic slip, supporting earlier earthquake-specific findings (Das & Henry, 2003), such as from 15 October 1979 M w 6.5 California (Doser & Kanamori, 1986;Hartzell & Heaton, 1986) and 11 March 2011 M w 9.1 Tohoku .  ; solid line depicts the aftershock regions from this study. Black lines depict active faults of Koehler et al. (2012). Magenta contours show the subduction interface for Slab 2.0 (Hayes et al., 2018), from 10 to 80 km in 10 km intervals. Thick red line is the trench (Krabbenhoeft et al., 2018). Event 869868, at the southwestern tip of the Kenai Peninsula, was ignored in drawing the aftershock region; event 868717, in western Kodiak Island, was not ignored (see Section 5). Considerations of coseismic slip, aftershock regions, and mainshock hypocenter motivate three main questions: 1. What is the relationship between coseismic slip and an aftershock region? 2. Is there a deficiency of aftershocks in the region of highest coseismic slip? 3. What is the relationship between the mainshock hypocenter and the regions of high or low coseismic slip? (Previous studies suggest there is none (Das & Henry, 2003;Henry & Das, 2001).)?
Potentially one could accumulate coseismic slip models for the 12 earthquakes-especially the most well-studied ones (like 2020 and 2021)-and compare them with the revised aftershock regions. To properly manage uncertainties, which are often not provided with coseismic slip models, this would be a major undertaking that is beyond the scope of our study.
To identify a specific challenge, consider a set of aftershocks (and mainshock) that we used to determine the aftershock region, such as in Figure 14. The hypothesis of Wetzler et al. (2017) suggests that any void of aftershocks in the interior of the aftershock region could possibly be associated with high coseismic slip. If a high coseismic slip occurred at the margin of the rupture area, it is possible that no aftershocks would occur and therefore the area would not be included in determining the aftershock region. This is a possible scenario for the near-trench easternmost portion of the 1957 aftershock region. In this area, known as the Unalaska gap  (see Section 6.7), there were no 1957 aftershocks to define the aftershock region, but there were tsunami observations from 1957 and paleotsunami deposits which suggest that shallow rupture has occurred (D. J. Nicolsky et al., 2016;.

Possibility of Near-Trench Slip on the Megathrust
Megathrust earthquake ruptures closer toward the trench tend to produce larger tsunamis (for a given earthquake magnitude). Kanamori (1972), focusing on the 1946 earthquake, defined "tsunami earthquakes" as events deficient in higher-frequency radiation and wrote: "The abnormal slow deformation of the tsunami earthquakes may be a manifestation of viscoelasticity of a weak zone beneath the inner margin of the trenches" (p. 346). Satake and Tanioka (1999) performed numerical modeling of the 1996 M w 7.5 Peru earthquake and concluded: "Most moment release of tsunami earthquakes occurs in a narrow region near the trench, and the concentrated slip is responsible for the large tsunami" (p. 467). The 2011 Tohoku earthquake revealed that tens of meters of lateral slip could occur near the trench (Fujiwara et al., 2011;Sato et al., 2011), countering the then-prevailing view that significant coseismic slip is not expected below the most updip portion of the megathrust (Hubbard et al., 2015;Pacheco et al., 1993), especially in accretionary subduction margins (Clift & Vannucchi, 2004).
From the aftershock regions (see Figures listed in column 2 of Table 5), we can get a basic idea of which ones are closest to the trench-and also which near-trench portions have not exhibited aftershock zones since 1938. The first item can be quantified, as shown in Table 5. For each aftershock region, we generate a set of uniformly distributed points within the region, then calculate the shortest distance to the trench of Bird (2003). This provides a distribution of distance-to-trench values, from which we calculate the minimum, mean, and maximum (Table 5 and Figure S24 in Supporting Information S1). We also calculate the distance d t below which 25% of the aftershock region resides. The 1946 aftershock region is exceptional for having the lowest d t (18 km).
The 1946 earthquake appears to have ruptured almost to the trench, based on the locations of the mainshock and aftershocks. It also produced one of the largest trans-Pacific tsunamis over the past century (J. M. Johnson & Satake, 1997;Okal & Hébert, 2007), with wave heights of 8 m in the harbor at Hilo, Hawaii (Macdonald et al., 1947). The revised aftershock regions for the other four earthquake-1964, 1938, 1957, and 1965-all Figure S1 in Supporting Information S1). (c) All events in this study that are within ±100 days of the mainshock. show a considerable distance between the aftershock region and the trench, in comparison with previous estimates ( Figure 1). Relevant questions regarding these near-trench areas are: 1. Is the interseismic plate motion aseismic, possibly including postseismic viscoelastic relaxation? 2. Did they rupture during the mainshock but not produce aftershocks? 3. If they have not recently ruptured, are they capable of rupturing in the near future as large earthquakes and producing large tsunamis?
Other data sets can help with these questions. For example, tsunami modeling may demonstrate the occurrence of high slip near the trench (D. J. Nicolsky et al., 2016), while plate coupling derived from GPS measurementsideally from seafloor geodetic sites-would inform strain accumulation estimates (Freymueller et al., 2008). It is possible that outer-rise events could be used to help inform the occurrence of near-trench slip (Sladen & Trevisan, 2018;Wetzler et al., 2017). Table 5 displays two additional calculations besides the proximity of each aftershock region to the trench: depth of the slip region (km) and average slip (m). If we assume that each aftershock region represents coseismic slip on the megathrust (Slab 2.0), then we can calculate the minimum, mean, and maximum depth of the slip region; an example for the 1946 aftershock region is shown in Figure S24 in Supporting Information S1. The mean depth is <21 km for 6 events (1938, 1946, 1964, 1965, 2003, and 2021) and >27 km for 5 events (1948, 1957, 1986, 1996, and 2020).

Depth of Aftershock Regions and Inferred Average Coseismic Slip
We also consider the average slip within each aftershock region that would be consistent with the published magnitude, under the assumption that the aftershock region is the region of coseismic slip. The calculation of average slip S is derived from seismic moment M 0 = μAS, where A is the area and μ is a shear modulus value for the crust. (We note that the magnitude estimates we use (Table 1) are not derived from the aftershock spatial extent but rather from modeling long-period seismic waves or tide gauge records of the tsunami). Average slip values range from 1 m (1986) to 12 m (1964). The 4 events with average slip >4 m all occurred at shallower depths (see "Slab distance" in Table 5) : 1946, 1964, 1965, and 2021. The 5 events occurring at greatest depths (1957, 1948, 1986, 1996, and 2020) all have average slip <4 m. It is possible that the higher-slip, shallower regions are associated with higher plate coupling than the down-dip regions (Elliott & Freymueller, 2020 Note. "Trench distance" is the shortest distance from a point inside the aftershock region to the trench. "Depth" is the depth (km below sea level) of a point interior to the aftershock region, projected onto the subduction interface. The variable d t is the distance from the trench containing 25% of the aftershock region area: P(d < d t ) = 0.25. The average slip S is directly calculated from the magnitude and area by using M 0 = μAS, where μ = 44.1 GPa is the shear modulus (PREM at 20 km depth) and seismic moment M 0 is calculated from M w . See Figure S24 in Supporting Information S1 for an example calculation.

Table 5
Calculations Derived From Our Aftershock Regions, the Location of the Trench (Bird, 2003), and the Subduction Interface (Hayes et al., 2018) (Section 6.2)

Aftershock Region of the 1964 Earthquake
The revised aftershock region for the 28 March 1964 M w 9.2 is different from the previously published aftershock regions of  and McCann et al. (1979) in several important places. (Note that the 1964 aftershock region of  is the same as that of McCann et al. (1979).) The most notable differences are in the eastern region, where our region extends about 120 km east of the region from McCann et al. (1979) (Figure 10, Table S1 in Supporting Information S1), into the Pamplona fault zone. The thrust mechanism of a 1964 aftershock (Stauder & Bollinger, 1966a) provides support for the nearby aftershocks occurring on a shallow-dipping extension of the megathrust. This is supported by additional studies of earthquakes in the Pamplona fault zone (Doser et al., 1997;Perez & Jacob, 1980).
The eastern extent of the aftershock region occurs in a region where the megathrust is known to be active during interseismic and seismic periods, based on current locking estimates Suito and Freymueller (2009);Elliott et al. (2013), coseismic displacement in 1964 (Plafker, 1969), and pre-1964 ruptures derived from paleoseismic data (tsunami deposits and coseismic vertical elevation changes) (Shennan et al., 2009).
The revised aftershock region provides a clearer understanding of the modern (1976-2022) thrust events occurring near the trench south of Kodiak ( Figure 9): they occur outside-not inside-the aftershock region of the 1964 earthquake.

Notable Foreshocks
Detailed studies of seismicity and geodetic data prior to the 11 March 2011 M w 9.1 Tohoku and 27 February 2010 M w 8.8 Maule earthquakes show that large earthquakes are preceded by transient signals lasting up to years Ruiz et al., 2014;Schurr et al., 2014;Yokota & Koketsu, 2015). Although we included foreshocks within 100 days in our initial event selection, we opted to exclude these when establishing the most likely events associated with the mainshock rupture. Here, we list some foreshock occurrences that could possibly be physically linked with the mainshocks. Figure 5 shows the events excluded from our analysis due to occurrences within 100 days prior to the mainshocks. Zooming in on this figure ( Figure S26a in Supporting Information S1), we find a region near Umnak Island (within the Fox Islands) that hosted 8 M > 6 earthquakes between 03 December 1956 and 09 January 1957, including six events on 02 January 1957. These events occurred two months prior to the mainshock and also hundreds of km east of the mainshock epicenter. Two of the events were analyzed by Stauder & Udias, 1963 and exhibit thrust mechanisms similar to other mechanisms of aftershocks for the 1957 earthquake. The region south of Umnak Island also produced earthquakes prior to the 10 June 1996 M w 7.9 earthquake ( Figure S26a in Supporting Information S1); the 1996 mainshock occurred hundreds of km west of Umnak.
Two smaller mainshocks within the 1957 aftershock region each experienced their own foreshocks. Two hours prior to the 07 May 1986 M w 7.9 earthquake there was a M w 6.2 foreshock on the southern (updip) boundary of the aftershock region. A day before the 10 June 1996 M w 7.9 earthquake there was a M w 6.5 foreshock in the western portion of the aftershock region; earlier events started with M w 5.  Figure S1 in Supporting Information S1). This foreshock event inspired the study of Stauder and Bollinger (1966a), who wrote: "After this discovery of the close parallelism between the [mechanisms of the] main shock and the large preshock, an effort was made to investigate as many as possible of the aftershocks" (p. 5287). We note that although the foreshock event occurred within 100 km of a large set of aftershocks south of Kodiak Island, it occurred outside the aftershock region and hundreds of km from the mainshock epicenter. (All foreshocks, including these two, were excluded from the final set of 324 events used to define the aftershock regions.)

Megathrust Earthquakes in the Early Instrumental Era (1897-1937)
There appears to have been no other M w ≥ 8 earthquake on the Aleutian-Alaska megathrust between the start of the instrumental era, in 1897, and the 10 November 1938 M w 8.3 earthquake. Examining the catalogs of Abe and Noguchi (1983) and ISC-GEM, we find 8 events with magnitudes above 7.5. Two are the 1899 Yakutat Bay events, which are not considered to have ruptured the megathrust (Plafker & Thatcher, 2008), one is the 17 August 1906 intraslab event (Okal, 2005) and one is a 1929 normal-fault event (Kanamori, 1972). The epicenter of the 14 June 1898 M s 7.6 earthquake is at the western extent of the 1965 region; its data were insufficient for the event to be relocated by Boyd and Lerner-Lam (1988). The 17 December 1929 M s 7.8 earthquake is also at the western extent of the 1965 region and also within the eastern extent of the 17 July 2017 M w 7.8, which had a similar magnitude . The 09 October 1900 M w 7.6-8.0 offshore Kodiak earthquake could have occurred on the megathrust, but other source regions are plausible (Tape et al., 2021). Therefore, from 1897 to 1937 the number of possible Aleutian-Alaska megathrust earthquakes with magnitudes 7.7-8.0 is between 0 and 3 (1898, 1900, and 1929).

Revisiting Sykes (1971)
With the benefit of an additional 50 years of earthquakes and also revised aftershock regions, we revisit . His three objectives were:

Summary
Aftershocks of major earthquakes provide valuable information for interpreting the slip areas of coseismic ruptures. This is particularly important for older earthquakes, for which data availability, quality, and station coverage may not be sufficient to constrain the details of a finite-source rupture model. We have relocated 553 earthquakes associated with 12 M w ≥ 7.7 earthquakes along the Aleutian-Alaska megathrust. We have used modern relocation methods that include probabilistic characterization of uncertainties, as well as systematic station errors. Much of our study documents our selection of events for the final catalog of 324 used to define the revised aftershock regions.
The aftershock regions of five earthquakes (1938, 1946, 1957, 1964, and 1965) span 3,548 km-nearly the full length of the megathrust-from the western end of 1965 to the eastern end of 1964 (Figure 1d). This was first documented in , and here we emphasize some new findings for these events based on the revised analysis.
1. The aftershock region for 1964 ( Figure 10) extends east into the Pamplona thrust region. This eastern limit of coseismic rupture is supported by coseismic vertical displacements and, indirectly, by the presence of interseismic coupling (Section 6.4). 2. The aftershock region for 1938 ( Figure 11) is almost twice as long as the one in McCann et al. (1979) (Table 4).
The revised region spans three other large events (1948, 2020, and 2021) and overlaps in the east with the western extent of the 1964 region. The region southwest of Unga Island experienced aftershocks in 1938 (Figure 11), 1948 (Figure 16c), and 2020 ( Figure 16a). South of this region, toward the trench, there is a gap of historical and modern earthquake occurrence (Figures S25d-S25f in Supporting Information S1). 3. The aftershock region for 1946 is remarkably close to the trench ( Figure 12) and implies that coseismic slip (and seafloor deformation) occurred on the megathrust close to the trench and generated the trans-Pacific tsunami (J. M. Okal & Hébert, 2007). The aftershock region is 1.3 times larger than that of McCann et al. (1979) and 2.4 times larger than that of   (Table S1 in Supporting Information S1). 4. The aftershock region for 1957 is 1,230 km long and is shifted west from the region of McCann et al. (1979).
Most importantly it is 0.76 times the area (Table 4), due to the absence of aftershocks along most of the updip section of the megathrust. 5. The aftershock region for 1965 is shorter (by 110 km) and 0.52 times the area of the region of McCann et al. (1979) (Tables 4 and S1 in Supporting Information S1). 6. The revised aftershock regions show that those for 1965 and 1957 are partially overlapping (Figure 15), as are those for 1938 and 1964 ( Figure 11) (Section 6.7). This has implications for interpreting the presence of structural barriers to earthquake ruptures (e.g., Briggs et al., 2014;Kaneko et al., 2010;Lay & Kanamori, 1981). For example, when determining scenario earthquakes (D. Nicolsky et al., 2013), it may be safest to ignore the notion of structural barriers.
Modern events provide the best opportunity for understanding the relationship between coseismic rupture areaderived from seismic and tsunami waveforms of the mainshock and geodetic observations-and the aftershock region. These events also provide an opportunity to understand the connection between interseismic locking and coseismic slip, something that was not possible prior to GPS measurements that became available in the 1990s. The 5 modern (post-1976) earthquakes we study (1986, 1996, 2003, 2020, and 2021) have all occurred as smaller events within the five historical aftershock regions: 2020 and 2021 within the 1938 region; 1986 and 1996 within the 1957 region; and 2003 within the 1965 region.
The aftershock regions of previous studies, beginning with , emphasized coverage of-and absence of-regions along the direction parallel to the trench; examples included gaps and abutting aftershock regions. Our revised aftershock regions enable interpretations in the trench-normal direction, based on quantities in Table 5. Our revised aftershock regions should be beneficial for understanding the connections between earthquake rupture zones and variations in plate coupling, volcanism, and geological structural variations within the upper and downgoing plates. The long list of studies in the introduction shows the value of such aftershock regions, and some of the findings of these studies may be impacted by the revised aftershock regions.
We have focused on only one part of the plate boundary system examined by  and later studies. The Queen Charlotte fault system (Brothers et al., 2020;St. Amand, 1957) would be a natural extension of this work. Future studies may identify some of the 324 Aleutian-Alaska events as not occurring on the subduction interface, either based on hypocenter or source mechanism. This reduction of events could lead to reductions in the size of the aftershock regions.  An overarching goal of seismologists and geophysicists is to understand the slip budget of the subduction system, that is, the fraction of slip occurring coseismically, postseismically, or aseismically for each portion of the Aleutian-Alaska megathrust. This information, combined with century-scale records of historical earthquakes and tsunamis, provides a foundation for interpreting the most likely sections of the megathrust to rupture in the next large earthquakes.

Data Availability Statement
• Our event selection is based on events within the ISC-GEM catalog, v. 9.0 (International Seismological Centre, 2022), which also contains a supplemental catalog. • Earthquake source mechanisms, represented by beach ball symbols, are plotted for the GCMT catalog (Dziewonski et al., 1981;Ekström et al., 2012) https://www.globalcmt.org/ and for the compilations by Stauder (Stauder & Udias, 1963;Stauder & Bollinger, 1966a;Stauder, 1968aStauder, , 1968b). • All earthquake relocations were performed with NonLinLoc (Lomax et al., 2000(Lomax et al., , 2014 (http://www.alomax. net/nlloc; https://github.com/alomax/NonLinLoc; last accessed February 2022). • Complete results from our NLL-relocated earthquake catalog of 553 events are available at Lomax (2022). This includes arrival time data, which were obtained from ISC, either from a public internet-based query or from files emailed to us. The 50 events with custom arrival times are identified in our master text file (Data Set S2) that lists our results. • A collection of digitized aftershock regions from previous studies (Figure 1) is available at Tape (2022). • Our own aftershock regions are available as digitized files (Data Set S1).