New criteria for selecting reliable Thellier-type paleointensity results from the 1960 Kilauea lava flows, Hawaii

Thellier-type paleointensity experiments associated with partial thermal remanent magnetization checks have been widely used to determine paleointensity values from volcanic and archaeological media. However, previous studies have revealed that a substantial portion of paleointensity results with positive checks for historical lava samples largely fails to predict known Earth magnetic field intensity values. To determine the fidelity of paleointensity values, conventional Thellier-type paleointensity experiments were performed on Kilauea lava flows that erupted in 1960. The positive partial thermal remanent magnetization checks for our results range from 30.28 ± 1.38 µT to 52.94 ± 1.89 µT. This strongly indicates that conventional paleointensity checks cannot guarantee the fidelity of paleointensity results, especially when the unblocking temperatures for the newly formed magnetic particles are higher than the treatment temperature. Therefore, in this study, to check for thermal alteration during heating, the temperature dependence of the hysteresis parameter measured at room temperature for the thermally treated samples was also measured. Our new results show that nearly all biased paleointensity values correspond to a ratio of the coercivity of remanence to the magnetic coercivity of > 3 and a chemical alteration index > ~ 10%, which indicates the strong effect of the domain state and thermal alteration on the fidelity of the paleointensity results. Our study provides feasible criteria to further improve the fidelity of paleointensity estimations.


Introduction
Variations in the Earth's magnetic field intensity at different time scales contain significant information concerning the growth of the Earth's deep interior (Macouin et al. 2004;Tarduno et al. 2007;Biggin et al. 2015) and the evolution of the geodynamo (Larson and Olson 1991;Olson et al. 2013). There are two different types of approaches used to trace temporal variations in the geomagnetic field intensity: relative paleointensity (RPI) determinations (Tauxe 1993) and absolute paleointensity determinations (Thellier 1938;Thellier and Thellier 1959). Sediments carry depositional or post-depositional remanent magnetization, and are excellent media to record the semicontinuous RPI (Valet et al. 2005;Yamazaki and Oda 2005;Channell et al. 2009;Ziegler et al. 2011). Meanwhile, volcanic rocks are suitable for high-accuracy spot readings of the ancient geomagnetic field intensity.
Modern absolute paleointensity determinations require multiple step heating with systematic alteration (or consistency) checks. The Thellier protocol (Thellier 1938) was initially proposed to compare the destruction of thermoremanent magnetization (TRM) and the acquisition of laboratory-induced partial thermoremanent magnetization (pTRM) at equal temperatures. Thellier-type double-heating techniques include multiple slight modifications, with each method having advantages and disadvantages (Thellier and Thellier 1959;Coe 1967;Aitken et al. 1988;Yu et al. 2004). The most commonly used technique is the so-called "Coe" protocol (Coe 1967) in which the specimen is first heated to an initial temperature (T i ) in a zero-field to determine the natural remanent magnetization (NRM) lost. To determine the pTRM gained, the specimen is reheated to T i under an in-field condition. Aitken et al. (1988) modified the Coe method (1967) by reversing the order of the double heating. The IZZI protocol alternates the Aitken method and the Coe method at every other step (Yu et al. 2004). Detailed reviews of various Thellier-type techniques have been provided by Valet (2003) and Biggin (2010).
Once a paleointensity determination is performed, the results are displayed on an Arai plot (Nagata et al. 1963) where the absolute value of the slope of the NRM remaining versus the pTRM gained indicates the ratio of the ancient to laboratory magnetic field intensities. Only stable single-domain (SD) particles follow linear Arai plots, reflecting identical spectra of the unblocking temperature (T ub ) and the blocking temperature (T b ). Multi-domain (MD) and pseudo-single-domain (PSD) grains tend to produce non-linear Arai plots (Dunlop and Özdemir 2001;Paterson et al. 2015). Inequivalent T ub over T b values are responsible for such non-linearities in Arai plots. The pTRM tail check was introduced to duplicate the previsouly acquired pTRM for non-uniformly magnetized TRM (Riisager et al. 2000;Riisager and Riisager 2001;Yu and Dunlop 2003). In addition to the physical origin, chemical contributions can alter the linearity of an Arai plot. For example, the growth of newly formed magnetic particles via chemical transformations (Yamamoto et al. 2003) can induce thermochemical remanent magnetization (TCRM).
The most convenient and easiest method to ensure high-fidelity paleointensity determinations is to check the quality of the paleointensity determination using historic rocks whose geomagnetic field intensities are readily known (e.g., using the International Geomagnetic Reference Field (IGRF), Thébault et al. 2015). Accordingly, the reliabilities of absolute paleointensity determinations have been tested using historic lava samples from Hawaii in the US (Tsunakawa and Shaw 1994;Chauvin et al. 2005;Oishi et al. 2005;Herrero-Bervera and Valet 2009;Morales et al. 2010;Cromwell et al. 2015), Italy (Calvo et al. 2002), Japan (Tsunakawa and Shaw 1994;Yu 2012), Spain (de Groot et al. 2015;Calvo et al. 2016), and the western US (Coe et al. 2004).
Of these various historic sites, the 1960 historic lava flows in Hawaii have drawn considerable attention because of their accessibility (Abokodair 1977;Tanaka and Kono, 1991;Tsunakawa and Shaw 1994;Tanaka et al. 1995;Valet and Herrero-Bervera 2000;Hill and Shaw 2000;Yamamoto et al. 2003;Chauvin et al. 2005;Hettero-Bervera and Valet 2009). However, relatively easier accessibility does not guarantee a successful duplication of the geomagnetic field intensity information. Böhnel et al. (2011) systematically summarized these paleointensity estimations and concluded that the paleointensity estimations depend on several factors including the lava cooling rate, the temperature treatment intervals, and experimental methods. In fact, the 1960 Hawaiian lava samples have often yielded paleointensity determinations biased toward higher/lower values by up to 10-20% or even more (Tanaka and Kono 1991;Tsunakawa and Shaw 1994;Tanaka et al. 1995;Valet and Herrero-Bervera 2000;Hill and Shaw 2000;Yamamoto et al. 2003;Herrero-Bervera and Valet 2009;Morales et al. 2010;Cromwell et al. 2015). For example, Tanaka and Kono (1991) investigated the 1960 lava flows on the Big Island of Hawaii and obtained a paleointensity that was approximately 50% higher than the expected value. Conversely, a recent microwave study using the IZZI + protocol produced paleointensity estimates consistent with the expected values (Grappone et al. 2019) and proposed that undetected chemical alterations are significant factors influencing the low paleointensity values found using conventional Thellier-type methods. The origin of heterogeneous or anomalous paleointensity outcomes likely lies in the influence of TCRM (Yamamoto et al. 2003), the presence of local magnetic anomalies (Morales et al. 2010), and/or the alteration and neo-formation of magnetic particles during repeated heating (Zhao et al. 2014). In addition to chemical alteration effects during paleointensity experiments, Fabian (2009) proposed that partial dissolution or low-temperature oxidization after TRM acquisition could seriously contaminate the original TRM but could not be distinguished from the true TRM in Thellier experiments. In addition, Cromwell et al. (2015) collected individual volcanic glass samples from the 1960 Hawaiian lava flows; these samples yielded nearly ideal estimation of the expected paleointensity, and had a much lower variance than published data from the same volcanic unit. This strongly indicates that domain effects should also be considered in paleointensity estimations.
The present study was designed to combine the domain state and the temperature dependence of magnetic hysteresis as a potential alteration check. We aim to determine the exact mechanisms for biased paleointensity results using positive pTRM checks. This is extremely important for paleointensity determinations of older rocks.

Samples
The Big Island of Hawaii, located on the Hawaiian hotspot, is currently volcanically active. The highest volcano on the island, Mauna Kea, is located to the north-west of the sampling site. Mauna Kea erupted from 300,000 to ~ 4000 years ago (Frey et al. 1990) (Fig. 1, upper red triangle). Mauna Loa is one of the largest volcanoes on the Big Island of Hawaii and is an active volcano, along with Kilauea. Kilauea is located at 19° 3′ N, 204° 3′ E (Fig. 1). On January 13, 1960, an eruption began along the east rift zone of Kilauea. The 1960 eruption lasted for approximately a month, and the lava covered an area of 10 km 2 (Ritcher et al. 1970;Rowland and Walker 1987). According to the IGRF model, the geomagnetic field intensity at the eruption site should have been 36.2 µT. Paleointensity determinations have yielded results ranging from 33.55 µT to 53.5 µT (Tanaka and Kono 1991;Tsunakawa and Shaw 1994;Tanaka et al. 1995;and Yamamoto et al. 2003). For example, Tanaka and Kono (1991) documented the presence of two segments in an Arai plot. They picked lower temperature segment because the pTRM checks failed at higher temperature ranges. As a consequence, their paleointensity determinations were overestimated (Tanaka and Kono 1991).
A total of 30 cylindrically shaped basaltic core samples were collected from the historic Kilauea flows using a portable engine drill in March 1998 (Yamamoto et al. 2003). The sampling site was a roadcut of exposed lava ~ 4 m thick and ~ 15 m wide. On the basis of the core positions, the samples were divided into four subsites, designated A1-A5, B1-B7, B8-B10, and C1-C4 (Yamamoto et al. 2003). The rock magnetic results (Yamamoto et al. 2003) indicate that the dominant magnetic carrier is the PSD low-titanium magnetite.

New rock magnetic measurements
Magnetic hysteresis measurements were performed using a variable field translation balance with a saturation field of 1 T. Values of the saturation magnetization (M s ), saturation remanence (M rs ), and magnetic coercivity (B c ) were calculated after removing the paramagnetic portion. The coercivity of remanence (B cr ) was determined using stepwise backfield demagnetization measurements. At least one chip was used for each of the 30 different basaltic cores. The magnetic hysteresis results are displayed on a Day diagram (Day et al. 1977) according to the criteria of Dunlop (2002) (Fig. 2a). All the samples fall within the pseudo-single-domain (PSD) range (Fig. 2a). The squareness versus coercivity plot shows a linear correlation with a mean slope value of approximately 0.01 mT −1 (Fig. 2b).
The first-order reversal curve (FORC) diagrams were also obtained (Smooth factor = 8), using an alternating gradient force magnetometer (MicroMag 2900). The results from a chip of the A1 sample display closing and elongated contours, reflecting the fine-grained nature of the sample (Fig. 3), as well as additional PSD features (Fig. 3).

Temperature-dependent rock magnetic properties
The paleointensity results can vary considerably as a result of several factors. Non-linear features in Arai plots can be caused either by MD effects or by thermal alteration during heating. To determine the exact reason for the non-linearity, sister samples were heated to the same paleointensity measurement heating steps (200, 300, 350, 400, 450, 480, 500, 520, 540, 560, 580, and 600 °C). Then, after heating, their hysteresis loops were measured at   Fig. 3: red indicates sample A1 (Fig. 3a); orange indicates sample B4 (Fig. 3b), magenta indicates sample B8 (Fig. 3c), green indicates sample C7 (Fig. 3d), blue indicates sample A11 (Fig. 3e), and purple indicates sample C2 (Fig. 3f ) room temperature for each step for each step to detect thermal effects (Haag et al. 1995;Henry et al. 2005;Smirnov and Tarduno 2003;Qin et al. 2011). The room temperature-normalized magnetic saturation remanence (Mrs room ) versus temperature is shown in Fig. 4. The M rs T /M rs room features can be divided into two groups from room temperature to 500 °C. The first group is relatively stable from 300 to 500 °C (Fig. 4a, e, f and g). The second group exhibits a dome-like feature ( Fig. 4b-d, h and i). For nature samples, the influences of the viscous remanent magnetization (VRM) usually reach 300 °C. Thus, to provide a more practical model, a new interval from 300 to 500 °C was chosen to recalculate the paleointensity (this is referred to as the restricted paleointensity) ( Table 1).
Paleointensity determinations were accepted when they passed the following selection criteria: (1) the demagnetization of NRM must be univectorial with the anchored maximum angular deviation (MAD) angles of < 4° (Kirschvink 1980); (2) all pTRM% checks must be within 5% (Paterson et al. 2014); (3) the angle between the principal fractions anchored to and free from the origin DANG (Tauxe and Staudigel 2004) should be < 4°; (4) at least six data points with more than 48% of the extrapolated NRM fraction (Coe et al. 1978) must be included in the regression analysis; and (5) the degree of scatter around the best-fit line normalized by the estimated slope (York 1969) should be < 5%. According to these selection criteria, a total of six samples had successful paleointensity estimations (Fig. 5; Table 1). The magnetization was removed at 560 °C (< 10%) with a few exceptions. The paleointensity results range from 30.28 ± 1.38 µT to 52.94 ± 1.89 µT. Several samples yielded biased values significantly higher than the expected values (IGRF model ≈ 36.1 μT). Such erroneous results with high estimated paleointensity have been attributed to the production of thermochemical chemical remanent magnetization (TCRM) during the heating processes (Yamamoto et al. 2003). A total of 6 out of 30 samples passed the stringent selection criteria (Table 1). The mean paleointensity of 43.40 ± 8.31 µT from 6 accepted samples disagrees with the predicted value (36.2 µT) by 19.9%.
For all six accepted samples, the temperature dependence of M rs indicates that the M rs values remain relatively stable from 300 to 500 °C (Fig. 4) and change significantly at higher temperature. Accordingly, we recalculated the paleointensity estimates using only the data points acquired from 300 to 500 °C (Table 1). Four samples yielded essentially identical intensity estimates within the uncertainty with or without the re-calculation (Fig. 6, Table 1). However, the paleointensity estimates from four samples (A1, B2, B5, and B7) were significantly lowered after the re-calculation using the 300-500 °C steps (Fig. 6, Table 1), which is more consistent with the  Table 1 Paleointensity results from the 1960 Hawaiian lava flows erupted from Kilauea ΔT is the temperature interval used in paleointensity estimation; n is the number of points used in paleointensity estimation; f, g, q, β and θ are NRM fraction, gap factor, quality factor, and the standard error of the slope and different angle of the selected NRM component from the origin on the orthogonal plot; Bg (μT) indicates the differences in the paleointensity estimations between the conventional IZZI approach and the more restricted values using the temperature interval from 300 to 500 °C; and * indicates that the value meets the condition (B cr /B c < 3 and CI < 10%) expected value although the experimental success rate is highly reduced.

Discussion
The studied lava flows erupted in 1960 from Kilauea on the Big Island of Hawaii. Altogether, 6 out of 30 samples yielded paleointensity results with positive pTRM checks using the IZZI protocol. The mean value of the paleointensity was ~ 43.4 µT, which is higher than the expected value for the study area (~ 36.1 µT). In addition, the paleointensity results range widely between 30.288 ± 1.38 and 52.94 ± 1.89 µT (Fig. 6).
Previous studies have already shown that positive pTRM checks do not guarantee the fidelity of paleointensity results for Thellier-type experiments because the unblocking temperatures of the newly formed magnetic particles could be higher than the treatment temperature (Coe 1967;Coe et al. 1978; Zhao et al. 2014). However, such a flaw can be overcome by examining the room temperature dependence of the magnetic parameters (Smirnov and Tarduno 2003;Qin et al. 2011). Upon thermal treatment, two types of mineral transformations occur. The first type of alteration can produce newly formed strongly magnetic minerals. Despite the distribution of the unblocking temperature, the enhanced concentration dependence of the magnetic parameters can sensitively detect the neo-formation of magnetic minerals.
The second type of thermal alteration does not produce new minerals, but can oxidize the primary magnetic minerals into a less magnetic state, e.g., from magnetite to maghemite, which could significantly decrease the saturation magnetization (also ARM and SIRM values). The first and second types of thermal alterations decrease and enhance the estimated paleointensity values, respectively (Zhao et al. 2014). Accordingly, Qin et al. (2011) adopted a new parameter called the chemical alteration index (CI), which can be defined, for example, by M rs_t °C /M rs_raw , where M rs_t °C and M rs_raw indicate the M s values for the samples after the t °C thermal treatment and the raw sample, respectively. The selection of the treatment temperature, t °C, should be determined according to the exact experiment results. For example, in this study, we selected 300 °C and 500 °C for CI to track the possible chemical changes and this temperature selection is more practical for natural samples. In addition to the conventional pTRM checks, CI values should be less than a threshold value (e.g., < 10%), to indicate weakly thermal alteration effects. Such a rationale has been applied in this study.
In addition to the thermal alteration effects, PSD/ MD particles can also seriously distort the Arai plot. To simultaneously consider these two effects, we constructed a new correlation between CI and B cr /B c (Fig. 7). Clearly, nearly all samples with higher paleointensity have B cr /B c > 3 and CI > 10%. This strongly indicates that both of these factors can distort the Arai plot, yielding biased paleointensity values. In contrast, samples with B cr /B c < 3 and CI < 10% yield an average paleointensity value of 39.70 ± 7.56 μT, which is closer to the expected value (Table 1). This strongly indicates that positive pTRM checks indeed cannot exclude a large amount of failed estimations, likely because of the mechanism proposed by Qin et al. (2011) and Zhao et al. (2014).
Our study provides practical new criteria to determine the fidelity of the paleointensity, especially for results with conventional positive checks. This raises a serious question for existing paleointensity results from the geological record. We strongly suggest that the domain state and thermal stability of samples be systematically incorporated into the paleointensity studies, except for all conventional criteria.

Conclusions
We presented the paleointensities of 1960 Kilauea lava flow samples using the temperature dependence of the hysteresis parameters according to modified Thellier experiments using the IZZI protocol. Only 6 out of 30 samples succeeded according to the positive pTRM checks. However, a large proportion of the "successful" paleointensity results Plot of the coercivity ratio versus the variation ratio of normalized saturation remanence at room temperature (M rs-room ) for the nine samples that passed the selection criteria failed to predict the expected value. This strongly indicates that positive pTRM checks cannot completely exclude failed results because conventional pTRM checks cannot detect thermal alterations when the unblocking temperature of the newly formed magnetic particles is higher than the treatment temperature. Therefore, we propose that a combination of the chemical index and the domain state parameters be used in future studies, especially for older samples, to improve the fidelity of the paleointensity results.