INTRODUCTION

The current configuration of the Earth's Magnetic Field (EMF) is well known from the data obtained by global magnetic observatories and satellite missions, describing the variation of both dipolar and non-dipolar components. The reconstruction of the past geomagnetic field is achieved using paleomagnetic records. Different types of variations can be distinguished regarding their magnitude, duration and the global or regional character. Apart from obtaining new and reliable paleomagnetic records from different geological times, it is also necessary to develop mathematical models of the behaviour of the ancient geomagnetic field, which help to understand the fine characteristics of the EMF.

The geomagnetic excursions and reversals are of particular interest in modern geomagnetism and paleomagnetism research. Polarity transitions are generally considered as an event of relatively short duration, usually spanning 10³ - 10⁴ years (e.g., ^{Merrill and McFadden 2003}). Excursions are defined in terms of very brief (<10³ years) deviation of virtual geomagnetic pole (VGP) positions from the geocentric axial dipole (GAD) that lies outside the range of secular variation for a particular population of VGP (^{Laj and Channell, 2007}).

The geomagnetic polarity time scale, obtained from numerous studies conducted around the world, is mostly based on marine and lacustrine sediments (^{Laj and Channell, 2007}). The paleomagnetic excursions may provide invaluable information on the behaviour of the geodynamo during the transitional state. It is particularly important for the Late Pleistocene because the refinements related with extinction events and human evolution during some critical periods (^{Goguitchaichvili et al., 2009}). Secular variation and field reversal rate are strongly influenced by the variations at the core-mantle boundary (^{Glatzmaier and Robert, 1997}). Volcanic rocks are considered as reliable paleomagnetic recorders because of the high stability of the thermoremanent magnetization that provides instantaneous record of the ancient EMF (e.g., ^{Prévot et al., 1985}), but the records are dis- continuous because of the sporadic character of volcanic eruptions.

The present study is aimed to contribute to the time-averaged field global database and geomagnetic polarity instability time scale for the last 5 Myr. For this purpose, we collected recently Ar-Ar dated lava flows associated to the Colima volcano.

GEOLOGICAL SETTING AND SAMPLING

The Colima volcanic complex (CVC) is a volcanic chain oriented N-S, and is composed by three andesitic stratovolcanoes: Cántaro, Nevado de Colima, and Colima. The CVC is located in the central part of Colima graben and belongs to the western portion of the Trans- Mexican Volcanic Belt (TMVB), one of the largest continental volcanic arcs on the American continent with more than 1000 km of length (^{Luhr and Carmichael, 1990}). This volcanic plateau, approximately 1,000-2,000 m high, roughly extends from the Pacific Ocean to the Gulf of Mexico (Figure 1).

Figure 1 Location of the Trans-Mexican Volcanic Belt (TMVB, green area). Tectonic environment and principal volcanic centers of the TMVB. The blue box rep- resents the area covered by the Colima volcanic complex (CVC) with its two main structures, Nevado de Colima and Colima volcano. Numbers at the trench refer to the age of the plate entering the subduction zone. Triangles indicate the main stratovolcanoes of the TMVB: The Navajas (LN); Sangagüey (SA); San Juan (SJ); Ceboruco (CB); Tequila (TQ); Volcano Colima (VC); Tancítaro (TA); Palo Huérfano (PH); La Joya (HO); El Zamorano (ZA); Nevado de Toluca (NT); Iztaccíhuatl (IZ); Popocatépetl (PO); La Malinche (LM); Cofre de Perote (CP); Pico de Orizaba (Por). The circles indicate volcanic calderas: La Primavera (LP); Los Azufres (Al); Amealco (AM); Huichapan (HU); Los Humeros (LH); MGVF is one of the principal volcanic fields: Michoacán-Guanajuato volcanic field; Chic is the Sierra de Chichinautzin; Pms Chiconquaco-Palma Sola. 

The Colima volcano has a large historical volcanological record due to its intense activity (^{Cortés et al., 2010}). The formation of the CVC main buildings started about 1.5 Ma with the building of the Cántaro volcano, and a relatively intense activity continued until approximately 1 Ma (^{James et al., 1986}). Later, the volcanic activity moved about 15 km to the south, with the formation of Nevado de Colima, which involves three periods of eruptive activity between 0.53 and 0.15 Ma (^{Robin et al., 1987}). The next important event took place at about 50 ka, 5 km to the south, with the formation of the Colima volcano, started by the building of the Paleofuego volcano (Robin et al., 1987), where several consecutive collapses occurred (^{Luhr and Prestegaard, 1988}; ^{Robin et al., 1987}; ^{Komorowski et al., 1997}, ^{Cortés et al., 2005}, ²⁰¹⁰).

Due the large and permanent activity of the Colima volcano until now, it has been the subject of numerous studies offering detailed historical records and some absolute dating using mainly K-Ar systematics (^{James et al., 1986}; ^{Robin et al., 1987}). The major effort for dating purpose is due to ^{Cortés (2015)}, who recently reported 30 new absolute Ar-Ar radiometric ages. In contrast, the eruptive history of the Cántaro and Nevado de Colima volcanoes is still relatively poorly constrained.

Our sampling strategy was based on the geological studies of the CVC by ^{Cortés (2015)} (Figure 2a). We sampled 21 out of the 30 sites (Figure 2b) reported in that work for the Colima and Nevado de Colima volcanoes, prioritizing fresh outcrops with no alteration and relatively easy access. No tectonic tilt correction was applied in this study since all studied lavas were found sub-horizontal. When possible, samples were distributed vertically and horizontally over several meters in order to avoid some local effects due of block tilting. On average, nine standard paleomagnetic cores were obtained from each cooling unit. They were obtained with a portable drill and oriented using both magnetic (Brunton) and solar compasses. In few cases, however only magnetic orientation was possible and local magnetic declination is considered as the correction factor.

Figure 2 Geological map of the Colima volcanic complex area showing sample site locations. a) All 30 available radiometric dates from the area (Córtes et al., 2015). b) Location of the 21 paleomagnetic sites for this study. 

IDENTIFICATION OF MAGNETIC CARRIERS

In order to identify the magnetic carriers responsible for the remanent magnetization and to obtain information about their paleo- magnetic stability, several rock-magnetic experiments were carried out. These experiments included:

1) The acquisition of magnetic susceptibility curves in low field as a function of temperature helps to determine the Curie temperatures of the main magnetic minerals by means of the differential method described in Tauxe (2010). Continuous (K-T) measurements in air were performed with a MS2 Bartington susceptibility bridge equipped with a furnace with a temperature range of 30 °C - 650 °C.

2) Magnetic hysteresis experiments. The hysteresis loops and as- sociated isothermal remanent magnetization (IRM) acquisition curves were measured using a variable field translation balance. Measurements were carried out on whole-rock powdered specimens, and in each case, first IRM acquisition and backfield curves were recorded first.

Typical results of rock-magnetic experiments are reported in Figure 3 (Sites Col 2, 3 and 7). In most of the cases the thermomagnetic curves reveals Curie temperatures of 560°C, the presence of Ti-poor titanomagnetite as the unique carrier of remanence, and indicate moderate degree of alteration due to heating. In some specimens, Tirich titanomagnetite seems to co-exist with the almost pure magnetite phase (sample 94C064A corresponding to site Col 7).

Figure 3 A summary of rock-magnetic experiments for the most representative samples from the Colima volcanic complex: Susceptibility vs. temperature curves (left), the red curve represents the heating curve and the blue is the cooling curve; hysteresis curves for an induced magnetic field (middle); and isothermal remanence acquisition curves (right) obtained with a variable field translation balance with the Day Plot (see also ^{Dunlop, 2002} a) to estimate the domain state (Single-Domain, SD; pseudo-single-domain, PSD; multidomain, MD; superparamagnetic, SP; MD-SD and SD-SP, mixture) of magnetic carriers of hysteresis parameters Mrs/Ms versus Bcr/Bc. a)Sample 94C034A from Col 3; b) Sample 94C064C from Col 7; c)Sample 94C015C from Col 2; and d) relation of the ratios of hysteresis parameters for the remaining samples. 

Corresponding hysteresis curves are symmetric yielding quite similar parameters, near to the origin, without evidence of wasp-waisted behaviour (^{Tauxe et al., 1996}), which probably reflect very restricted ranges of the opaque mineral coercivities. When judging the ratios obtained from the hysteresis curves, it seems that the samples have a pseudo simple domain PSD in the Day plot provided by ^{Dunlop (2002}) (Figure 3d).

Isothermal remanence acquisition curves are sensitive to the magnetic mineralogy, concentration and grain size properties. Almost all samples are saturated at about 300 mT applied magnetic field, which indicate the presence of a ferromagnetic phase with moderate coercivity as may be expected from magnetite and titanomagnetite grains (^{Tauxe, 2010}).

REMANENCE PROPERTIES

Remanent magnetization was measured using a JR-6 spinner magnetometer. At the initial stage, three specimens (belonging to different cores) per site were selected for detailed thermal and alternating field (AF) treatments in order to choose the most suitable demagnetization method. An ASC TD-48 furnace was used during the thermal treatment, while a Molspin AF-demagnetizer allowed sample demagnetization to 5mT up to 95 mT. The components of the remanence for each specimen and the site-mean paleomagnetic directions were determined by the method of the principal component analysis (^{Kirschvink, 1980}) and Fisher statistics (^{Fisher, 1953}).

In most cases, a stable single component was detected (Figure 4b, 4d and 4e; sites Col1 sample 94C001A, site Col4 sample 94C035A and site Col8 sample 94C062A), accompanied by a negligible viscous overprint. In a few cases, however, the presence of secondary compo- nents (Figure 4a, 4c and 4f that correspond to Col3 sample 94C026A, Col2 sample 94C011B, and Col10 sample 94C098A, respectively) are observed probably due to viscous magnetic overprint and easily removed. Around 50 samples were thermally demagnetized. However, AF demagnetization was found to be a more efficient cleaning method as may be evidenced for samples 94C011A and 94C011B. It should be also noted that alternative field treatment has a little limitation because in a few samples (example, specimen 94C058A, site Col 6) no complete demagnetization is obtained applying maximum available peak field of 95 mT. However, the determination of characteristic remanence components may be achieved unambiguously for these samples using the principal component analysis (^{Kirschvink, 1980}). Site-mean paleo- directions were determined for all sites (Table 1). These directions are quite precisely determined since in all cases the values of α95 are less than 10° which is common for volcanic outcrops.

Figure 4 Examples of orthogonal vector plots, showing thermal demagnetization since ARN0 with steps of 50°C until 540° (label with an initial "P" example P0540) for the. alternating field demagnetization since ARN0 with steps of demagnetization 5mT until 100mT (label with an initial "H" example H100). a) Sample 94C026A corresponding to the site Col3. b) Sample 94C001A, site Col1. c) Sample 94C011B, site Col2. d) Sample 94C035A, site Col4, including two orthogonal vector plots of events belonging to transitional geomagnetic regime. e) Sample 94C062A, site Col8. f) Sample 94C098A, site Col10. 

MAIN RESULTS AND DISCUSSION

Nineteen lava flows yielded a normal polarity magnetization, while two sites gave clearly defined transitional paleodirections (Table 1). Both transitional lavas were radiometrically dated. The paleodirections from site Col8 (Table 1), dated as 30±12 ka, correspond to the transitional geomagnetic regime. Tentatively, it may correspond to the Mono Lake (^{Benson et al. 2003}, ^{Negrini et al., 1984}) or Laschamp (Denham and Cox, 1960, ^{Liddicoat and Coe, 1979}) excursion, accord- ing to the available Ar-Ar radiometric ages. The Laschamp excursion was the first reported geomagnetic excursion, and is certainly the best known event in the Brunhes Chron (Chaîne des Puys, Massif Central, France; ^{Bonhommet and Babkine, 1967}). The mean paleodirection of Col8 is based on only four out of nine samples demagnetized. However, the directions are grouped yielding VGP latitude of about 42°, strongly deviated from the GAD (geomagnetic axial dipole) directions. These transitional directions can be correlated with the Mono Lake or the Laschamp events, which are dated at 28 ka and 40-45 ka, respectively. The Mono Lake event may be considered the best candidate because it is usually found in North America (^{Negrini et al., 2014}; ^{Benson et al., 2003}). Site Col10 shows as well defined transitional magnetic polarity and may be related to the Calabrian Ridge I and Portuguese Margin events, both found in marine sediments: in the Ionian sea (^{Langereis et al., 1997}) for the Calabrian Ridge I and in the north-east Atlantic Ocean for the Portuguese Margin (^{Thouveny et al., 2004}; ^{Carcaillet et al., 2004}). Site Col10, dated as 300 ± 90 ka, yields VGP latitude of -38° pointing to the intermediate geomagnetic regime; this can be correlated to the Portuguese margin, located in the North Atlantic Ocean, with has an age about 290 ka documented by ^{Thouveny et al. (2004)} and ^{Carcaillet et al. (2004)} using marine sediments, and correla ted to The Calabrian Ridge I (^{Langereis et al., 1997}), an excursion event located in the temporal window around 315-325 ka.

Table 1 Results of flow mean paleodirections for Colima volcanic complex lavas. Site: Name of the sample sites as described in ^{Cortés (2015)}; Coord: Geographical coordinates of the sampled sites; Age: Ar-Ar ages of the flows in ka; n/N: Number of specimens used (n) from the total of the specimens sampled (N); Directions: Flow-mean inclination (In) and declinations (Dc); k and α95:precision of parameters of Fisher statistics; VGPs: Virtual Geomagnetic Pole positions; Pol: magnetic polarity, where N is normal and T is transitional; K is the dispersion of the VGP, and A95 is the quality factor of the VGP. 

Traditionally, sites with low VGP latitudes (Figure 5) are removed from time average field (TAF) and paleosecular variation (PSV) studies for recent times (<5 Myr) (^{Johnson et al., 2008}). Generally speaking, the scatter of the VGP obtained should characterize the paleosecular variation (PSV) of the geomagnetic field for the given latitude and age (^{Cox, 1969}), but at the same time the dispersion may be biased by the transitional data which represents an excursions of the geomagnetic field. Thus, these data must be rejected from the population before to attempt to characterize the PSV. Transitional data commonly are not taken into account by the use of a conventional cut-off angle of 45° and 60° (^{Watkins, 1973}) in order to separate the paleosecular variation and transitional regimen (^{Johnson et al., 2008}).

Figure 5 Mean directions, inclinations, declinations and paleolatitudes of the 21 volcanic flows from the Colima volcanic complex, against age underlying the duration transitional events in these periods of time. 

The mean paleomagnetic direction obtained in this study, rejecting two transitional sites, is: Dm = 1.2°, Im = 38.2° (N = 19, α₉₅ = 5.5°, k = 39.1) (Figure 6a). The corresponding paleomagnetic pole (Figure 6b) position is λ_p = 270.1°, φ_p = 87.6° (A₉₅ = 5.3°, K = 35.9). The obtained direction is very close to the expected direction D_BC = 3.8°, I_BC = 38.1°, and D_T = 1.6°, I_T = 35.1° obtained for the last 5 Myr (Pliocene) and 10 Myr (part of the Miocene), as derived from the available reference poles for the North American craton (^{Besse and Courtillot, 2002} and ^{Torsvik et al., 2012}, respectively), and no rotation or flattening is present.

Figure 6 a) Equal area projection of the mean paleomagnetic directions obtained in this study; b) projection of the VGPs together with reference paleomagnetic poles recalculated from the stable North America (^{Besse and Courtillot, 2002}; ^{Torsvik et al., 2012}). 

An important issue discussed during the last decade in paleomagnetism is the relationship between the latitude and the VGPs scatter; this implies that at higher latitudes the scatter increases (^{McElhinny and Mc Fadden 1997}; ^{Johnson et al., 2008}; ^{Linder and Gilder 2012}). The idea of the latitude dependence of VGPs scatter (^{Cox and Doell, 1960}) depends critically on a set of data from about 20° of latitude as argued by ^{Johnson et al. (2008)}. The angular dispersion:

is the formula used to estimate the paleosecular variation, where:

(^{Cox, 1969}) is the total angular dispersion; N is the number of the sites used, δ _i is the angular distance of the i ^th virtual geomagnetic pole (VGP) from the axial dipole, SW is the within-site dispersion, and n- is the average number of samples per site. The dispersion obtained in this study is S _F = 13.0°, the lower confidence limit S _L = 12.7°, and the upper confidence limit SU = 15.2° using the calculation method of Cox (1969). The values agree relatively well (within the uncertainty) with the model G of ^{McElhinny and McFadden (1997)} for the last 5 Myr (Figure 7). However, our data agree almost completely with the TK03 curve proposed by ^{Tauxe and Kent (2013)} but are different to the model G of ^{Johnson et al. (2008)}, mainly due to the different selection criteria using a cut-off value dependence on dispersion parameter k. Moreover, in the compilation of ^{Johnson et al. (2008)}, a single reference from Mexico is used (^{Mejia et al., 2005}) with only 16 mean directions reported mainly from the Matuyama chron. The dispersion parameters obtained in the present investigation are compatible to the sites of similar latitude, such as Réunion and the South Pacific but a little bit higher than those from Hawaii for Bruhnes (0-0.78 Ma, ^{Lawrence et al., 2006}). In Mexico, Mejia et al. (2005) studied the secular variation for the Holocene reporting a mean direction of Dm=358.8°, Im=31.6°, α ₉₅=2.0° and k=29° and a scatter of 12.7 ° and with lower and upper 95% confidence limits of 11.9° and 14.1°. ^{Conte-Fasano et al. (2006)} conducted similar analysis of the Michoacán Guanajuato volcanic field reporting almost similar paleodirections: Dm=357.9°, Im=28.4, α ₉₅=7.3°, k=21, while the VGP scatter was estimated as 15.4° with lower and upper 95% confidence limits of 19.6° and 12.7°. These values are slightly higher than those obtained in the present study and even to the values reported in the model G of ^{Johnson et al. (2008)}. The value of the dispersion of VGPs obtained in this study is still consistent with the values expected from the models G and Tk03, 13.6° and 12.9°, respectively, however it is low compared with the model G of ^{Johnson et al. (2008)} with a value of 14.5.

Figure 7 Scatter of the VGP as function of the latitude for the last 5 Myr, compared with the model G of ^{McElhinny and McFadden (1997)} (blue line labelled McFadden, 1997), ^{Johnson et al. (2008)} (black line labelled Johnson, 2008) and the Tk03 model proposed by ^{Tauxe and Kent (2013)} (green line). 

CONCLUDING REMARKS

The flow-mean directions obtained in this study may be considered to be of primary origin (characteristic remanence). A full battery of rock-magnetic experiments shows that the magnetization is carried in most cases by Ti-poor titanomagnetite, probably resulting from oxy- exsolution of original titanomagnetite during the initial flow cooling. In addition, relatively high unblocking temperature spectra and moderate to high coercivities point to pseudo-single domain magnetic structure grains as responsible for remanence.

The paleodirections are rather precisely determined for all 21 analyzed sites, yielding relatively low within site dispersion. All sites point to normal polarity magnetizations as should be expected for the cooling units erupted during the Bruhnes chron. Two sites how- ever yielded clearly defined intermediary paleodirections that may be correlated to Mono Lake or Calabrian Ridge I short geomagnetic excursions respectively.

The mean paleomagnetic direction obtained in this study, rejecting two transitional sites, is: Dm = 1.2°, Im = 38.2° (N = 19, α ₉₅ = 5.5°, k = 39.1°). The corresponding paleomagnetic pole (Figure 7) position is λ_p = 270.1°, φ_p = 87.6° (A ₉₅ = 5.3°, K = 35.9°). These directions are practically undistinguishable (Figure 6a and 6b) from both the spin axis and the expected Plio-Quaternary paleodirections, as derived from reference poles for the North American craton (^{Besse and Courtillot, 2002}; ^{Torsvik et al., 2012}). This may indicate that no major regional tectonic rotation occurred in the area since about Pleistocene.

The dispersion parameters obtained in present investigation are compatible to the sites of similar latitude like Réunion and the South Pacific but a little bit higher than Hawaii for the last 1 Myr and agree well to the secular variation model of ^{Tauxe and Kent (2013)}.