ARTICLE IN PRESS Acta Astronautica 66 (2010) 1008–1016
Contents lists available at ScienceDirect
Acta Astronautica journal homepage: www.elsevier.com/locate/actaastro
Neural network and wavelets in prediction of cosmic ray variability: The North Africa as study case Neı¨la Zarrouk , Raouf Bennaceur ´ des Sciences de Tunis, Tunisia Laboratoire de Physique de la Matie re Condense´e, Faculte
a r t i c l e i n f o
abstract
Article history: Received 15 December 2008 Received in revised form 23 July 2009 Accepted 21 September 2009 Available online 17 October 2009
Since the Earth is permanently bombarded with energetic cosmic rays particles, cosmic ray flux has been monitored by ground based neutron monitors for decades. In this work an attempt is made to investigate the decomposition and reconstructions provided by Morlet wavelet technique, using data series of cosmic rays variabilities, then to constitute from this wavelet analysis an input data base for the neural network system with which we can then predict decomposition coefficients and all related parameters for other points. Thus the latter are used for the recomposition step in which the plots and curves describing the relative cosmic rays intensities are obtained in any points on the earth in which we do not have any information about cosmic rays intensities. Although neural network associated with wavelets are not frequently used for cosmic rays time series, they seems very suitable and are a good choice to obtain these results. In fact we have succeeded to derive a very useful tool to obtain the decomposition coefficients, the main periods for each point on the Earth and on another hand we have now a kind of virtual NM for these locations like North Africa countries, Maroc, Algeria, Tunisia, Libya and Cairo. We have found the aspect of very known 11-years cycle: T1, we have also revealed the variation type of T2 and especially T3 cycles which seem to be induced by particular Earth’s phenomena. & 2009 Elsevier Ltd. All rights reserved.
Keywords: Cosmic rays Neural networks Wavelets
1. Introduction Galactic cosmic rays come mainly from outside the solar system. The flux of GCRs that reach the surface of the Earth is influenced by solar activity, latitude, altitude, diurnal cycle. During high solar activity, emissions of matter and electromagnetic fields from the sun namely solar wind increase making it difficult for GCRs to penetrate the inner solar system and then reach the Earth. The GCR intensity is low when the solar activity is high and vice-versa constituting an approximately 11 year periodicity [8].
Corresponding author. Tel.: þ216 23066470; fax: þ216 71885073.
E-mail address:
[email protected] (N. Zarrouk). 0094-5765/$ - see front matter & 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.actaastro.2009.09.023
Beginning from the cosmic ray variation curves given by 10 different neutron monitors stations [1–7], we analyzed the variation functions and decomposed them in Morlet wavelets for a first phase and we have proceeded for reconstruction to zoom and predict periods details or singularities hidden behind the original curves describing cosmic ray variation in time space. We have used the real time data series of cosmic rays variabilities of 10 stations [1–7] from the network of neutrons monitors stations for the Morlet wavelets analysis, then knowing the characteristics of the different NM stations we have built the training inputs base of neural network system. Two kinds of training were performed with this inputs base. On the first side the outputs are the Morlet coefficients decomposition and the maxima of Morlet coefficients decomposition which can be then obtained for any position of a virtual NM on
ARTICLE IN PRESS N. Zarrouk, R. Bennaceur / Acta Astronautica 66 (2010) 1008–1016
the Earth, we have taken Tunis, Cairo, Alger, Rabat and Libya as examples of studies especially that the galactic cosmic rays variations were not investigated for these locations. In the other side we have then used the decomposition coefficients for Morlet reconstructions of cosmic rays variabilities curves in each selected point on the Earth.
2. Theoretical methods and data 2.1. Different stations In July 1997 the Moscow NM station was the first one, in the world presenting real time data to the internet, and then the number of other stations increases operating in various latitudes around the world. There are now 25 stations providing real or quasi real time data (in digital and/or graphical forms Mavromichalaki et al. [9]) elaborated automatically to ensure compatibility with other stations. The use of all stations as a unified multidirectional detector repair the accuracy of the measurements and made them higher (o0.1% for hourly data). This system has been designed with the capability to support a large number of stations. It’s worth noting that the collection system is able to provide reliable data, based on the issue that there are independent programs collecting simultaneously data from different stations in different ways. Nowadays there are 21 stations from which the described system collects data. The early detection of earth directed solar energetic particles (SEP) event by NMs offers a very good chance of preventive prognosis of dangerous particles flux and can alert with a very low probability of false alarm.
2.2. Neural network and wavelets Neural networks are used to model complex relationships between inputs and outputs or for a specific application such as pattern recognition or data classification, through a learning process. Neural network models can be used to infer a function from observations and also to use it this is then the utility of artificial neural network. Neural network is also particularly useful in applications where the complexity of the data or task makes the design of such a function by hand impractical [10]. Within neural cell we perform an application from our input data to output data in the form of matrix relation; this relation can reproduce hidden non linear relationship. Thus we have used neural network system to model relationships between characteristics of neutron monitors on a side and relative cosmic rays intensities or in another term decomposition coefficients for desired wavelet reconstructions on the other side. In this study multilayer perceptrons (MLP) are used. The architecture of MLP consisted of four layers, the input, both hidden and the output layers. The activity of the input units represents the raw information that is fed into the network, while nodes in hidden and output layers processed information. The behaviour of the output units
1009
depends on the activity of the hidden units and the weights between the hidden and output units. The nodes in hidden layers had no prescribed initial values and helped to allow complex relationships between the input and output nodes to evolve. The activity of each hidden unit is determined by the activities of the input units and the weights on the connections between the input and the hidden units. Information was transported from the input layer to the output layer by calculating the sum on each node, which was derived from combining all the nodes in the previous layer. During the training, the corresponding known outputs of the system were held in the output nodes to compare with the results produced by the network. Neural networks draw their power from modelling of their capacity to collect the high level dependences, implying several variables at the same time. A problem of over-training appears when the neural network system learns integrally, it will give bad results when a little different data are presented to him. They exist methods to optimize the phase of training so that the phenomenon of over or under training disappears, of which the technique of the ‘‘early stopping’’ and the regularization. The capacity of the neural network is controlled basically by the number of hidden units nh, then if it is too small the network has not enough parameters and cannot collect all the dependences which are used for modelling and predicting the values of the observed process [15]. If one chooses reciprocally, a too large value for nh, the number of parameters of the model increases and it becomes possible, during the phase of optimization of the parameters, to model certain relations which are only the fruit of statistical fluctuations characteristic of training unit used rather than of the fundamental relations of dependence between the variables. Being able to model any function if the number of hidden units is adequate. Thus the neural networks are universal approximators. To make sure that the neural network is satisfied by the fundamental relations of dependence, one adds in addition to training unit another unit called ‘‘validation unit’’ so that at the end of each training, one measures not only the error of training but also the error of validation. Once the optimization phase of the parameters is completed this error of validation is calculated. After trying different training models each one with a different number of hidden units, one can compare the training and validation errors. The training error decreases as the number of hidden units increases. The validation error, as for it, is high when the number of hidden units is weak, decreases with the increase in the number of hidden units, reached a minimum for a certain optimal number of hidden units, and then increases again when the number of units becomes too large. We thus adjusted the optimal number of hidden neurons to realize this compromise between the validation error and training error, what enabled us to model our functions with a good convergence of training phase. The wavelet theory involves representing general functions in terms of simple, fixed building blocks at
ARTICLE IN PRESS 1010
N. Zarrouk, R. Bennaceur / Acta Astronautica 66 (2010) 1008–1016
different scales and positions. We have used translations and dilations of one fixed function for wavelet expansion. Sophisticated wavelets are more powerful in revealing hidden detailed structures for example Morlet wavelet has been used to examine the processes, models and structures behind the variability of cosmic rays intensities [11]. Using the same wavelet expansion we present in this work a mathematical zoom to discover the hidden structures and to extrapolate unknown aspects in cosmic variation, especially for points on the Earth where informations concerning cosmic rays intensities variations are incomplete or absolutely absent. The Morlet wavelet is defined as a complex sine wave, localized with a Gaussian. The frequency domain representation is a single symmetric Gaussian peak. The frequency localization is very good. This wavelet has the advantage of incorporating a wave of a certain period, as well as being finite in extent.
3.1. Preparation of training samples: decomposition in Morlet wavelets The method implemented employs a multi level decomposition scheme via the Morlet wavelet transform, and interactive weighted recomposition. Different reconstruction strategies are discussed. The frequency domain transform of a real wavelet is symmetric around frequency 0 and contains two peaks. The nature of an oscillation with the CWT spectrum will vary greatly with whether the wavelet is real or complex. The complex wavelet will evidence a constant power across the time duration of the oscillation. Otherwise, a real wavelet produces power only at those times where the oscillation is at an extreme or where a sharp discontinuity occurs. When decomposing a non linear time series into time-frequency space, wavelet analysis is a useful tool both to find the dominant mode of variation and also to study how it varies with time [12]. The wavelet transform of a function y(t) uses spatially localized functions given by the coefficients of decomposition real and imaginary components of the Eq. (1) which were calculated and plotted versus scale and translation parameters j and k Z
tf ti
wIj;k ¼ j1=2
Z
R yðtÞgj;k ðtÞ dt; tf
ti
I yðtÞgj;k ðtÞ dt;
R gj;k ðtÞ ¼ g R
Stations
Study periods
Latitude (1 dl)
Longitude (1 dl)
Altitude (m)
Moscow Jungfraujoch Irkutsk Oulu Lomnicky´ Climax Sanae Potchefstroom Tsumeb Hermanus
1958–2007 1992–2007 1958–2007 1964–2007 1982–2008 1953–2007 1977–2002 1972–2002 1977–2002 1973–2002
55.47N 46.55N 52.28N 65.05N 49.20N 39.37N 70.19S 26.41S 19.12S 34.25S
37.32E 7.98E 104.02E 25.47E 20.22E 106.18W 02. 21W 27. 06E 17.35E 19.13E
200 3570 475 15 2634 3400 52 1351 1240 26
by (Eq. (2)). t2 gðtÞ ¼ exp io0 t 2
ð2Þ
We have used measured, limited and discrete time series, 49 years period of study for example for Moscow NM, thus we discretize expression of Eq. (1) as follows for Moscow NM:
3. Results and discussions
wRj;k ¼ j1=2
Table 1 List of neutron monitors investigated in calculations.
tk j ð1Þ
where ½ti ; tf is the study time interval. j is the scale dilation, compressing and stretching of the wavelet g used to change the scale, j determines the characteristic wavelength; k is the translation parameter, the shifting of g used to slide in time [13], and g* the complex conjugate of g. Morlet wavelet is a complex sine wave multiplied by a Gaussian envelope and given
wRj;k ¼ j1=2
49 X
R yðtl Þgj;k ðtl Þ
ð3Þ
l¼1
As we have operated in our previous work [14] for cosmic rays variation versus time given by the curve relative to Moscow neutron monitor. We have investigated data corresponding to yearly cosmic ray variations measured by 10 neutron monitors stations: Moscow, Climax, Oulu,Jungfraujoch, Irkutsk, Lomnicky´, Sanae, Potchefstroom, Tsumeb and Hermanus among the twenty five known stations. These data constitute the training samples we have reproduced in main details all curves of yearly cosmic ray variations relative to different stations then we have decomposed the different cosmic ray variations functions in Morlet wavelets. The list of neutron monitors used for calculations is brought in Table 1. 3.2. Training Training on main periods: In a first step and once we have decomposed the different cosmic ray variations relative to 10 different stations in Morlet wavelets, we have also derived the main periods studying and plotting the decomposition coefficients versus time, then we have notice the main three periods T1, T2, and T3 corresponding to three first maxima coefficients for the 10 stations. We present the periods for each station as outputs to neural network system, the inputs are geographic coordinates of stations. Training on decomposition coefficients: Once the tables of imaginary and real components of Morlet decomposition coefficients corresponding to each station are obtained, we present a certain number of coefficients necessary for reconstruction as outputs to neural network system. The matrix of real coefficients is reduced to a
ARTICLE IN PRESS N. Zarrouk, R. Bennaceur / Acta Astronautica 66 (2010) 1008–1016
Fig. 1. The main periods T1,T2, and T3 variations versus latitude longitude and for an altitude of 26 m.
Fig. 2. CR Tsumeb variation corresponding for the interval of minimum for T3 and CR Moscow variations nearest to the maximum of T3.
1011
ARTICLE IN PRESS 1012
N. Zarrouk, R. Bennaceur / Acta Astronautica 66 (2010) 1008–1016
vector for a single value of scaling parameter; the inputs are the geographic coordinates of corresponding neutron monitor stations. After the derivation of values on the output nodes, the latter were compared with the desired values and a back propagation algorithm was applied to adjust the weights to decrease the difference between the actual and desired predictions. The process was repeated iteratively using all cases in the training set until it met the least mean square error (MSE) between the target and actual output values. We have reached a precision of 109, this work present the training phase of our neural network study.
Table 2 Neural network main periods for North Africa cases. Main periods in years
Rabat
Alger
Tunis
Tripoli
Cairo
T1 T2 T3
18.01 24.50 82.55
14.04 21.39 58.56
12.42 20.90 46.33
12.19 21.12 44.29
12.13 20.05 58.35
We have used two training stages. The training test of which samples are all examples of stations except one example on which we try the credibility of test loading the results of training and then, comparing the reconstructed signal to observed signal of corresponding station. On the other side the main training in which the complete stations unit was used as samples unit, the results of this training are used for our aim of extrapolating numerically neutrons monitors’ stations network. 3.3. Prediction of main periods Once the system trained on main periods and in the same training interval for latitude and longitude, we have interpolated and plotted in 3D dimensions space for each couple of (lat, long) T1, T2, and T3 versus latitude, longitude and for a fixed altitude of 26 m which may describe any station or position on the Earth, the results of interpolation are as shown in Fig. 1. The main known period of 11–12 years modulating the curves of relative cosmic rays intensities is illustrated and found here. In fact the values of T1 calculated with neural network system are ranging from 9.52 to 15.94 when
Fig. 3. T1,T2, and T3 variations with altitude for a longitude of 0.0348 rd and a group of latitude values varying from 0 to p/2 rd.
ARTICLE IN PRESS N. Zarrouk, R. Bennaceur / Acta Astronautica 66 (2010) 1008–1016
varying latitude and longitude. The less important period T2 found with neural network system and corresponding then to a lower decomposition coefficient according to periods training is found fluctuating between 11.58 and 38.56. For T3 the period which must correspond to the
1013
less frequent cycle seems to oscillate here from 40.82 to 80.64 years and we notice a certain singularity occurring in the interval of longitude between 121 and 231 in which we can intercalate the Tsumeb NM station longitude then for which we have found indeed the lowest value of T3, 33
Fig. 4. Morlet reconstructions tests on the right compared to measured signals for cosmic rays variations on the left.
ARTICLE IN PRESS 1014
N. Zarrouk, R. Bennaceur / Acta Astronautica 66 (2010) 1008–1016
years with respect to other stations of training. This singularity occurs as a fluctuation for T1 but it is deeply seen in T3, it is probably related to Earth’s atmosphere phenomena. Indeed in the context of finding incident cosmophysical periodicities that may modulate terrestrial phenomena, little attention has been given to the large scale climatic phenomena: the Atlantic multidecadal oscillation (AMO). Common periodicities were analysed [16,17] between phenomena presumably associated to hurricanes : the Atlantic multidecadal oscillation (AMO) and the sea surface temperature (SST) versus cosmic rays, and on the other hand cosmic rays versus Atlantic hurricanes. A common period of 3072 years was found for total Cyclonal Energy, the total number of Tropical Storms landing in the Atlantic coast of Mexico and others. Thus these terrestrial phenomena indexes modulate or at least influence at the first order the less frequent period due to cosmic rays modulation. We show then a competitive mechanism measuring the relevance of cosmophysical phenomena with respect to terrestrial sources of affectation. It is worth noting that we have used for the training the period T3 found by Morlet decomposition, the study and training periods of different NMs stations were different thus we expect that T3 do not appears for all South Africa cases of applications or appears with the presence of certain singularity. The maximum of T3 is reached for values of longitude between 401 and 571, the longitude of Moscow NM of 371 is the nearest. From the curves of cosmic ray variations of Moscow and Tsumeb NMs we can see the difference of variations values for the cosmic rays intensities (Fig. 2). We have on another hand studied the variation in 2D dimensional space of main periods T1, T2, T3 separately with altitude, this is shown in Fig. 3. The most dominant period T1 decreases with altitude and also with latitude parameter for a fixed longitude and for all the latitudes varying from 0 to p/2 rd keeping also the same shape of variation. T1 is ranging from 8.08 to 15.17 years. The decrease of the most dominant period T1 with altitude is
expected. Indeed this known cycle of 11 years mean value is mostly related to extragalactic component of cosmic rays which increase in intensity with altitude and is the most energetic and important near the top of the atmosphere. This is the atmospheric shielding provided by the earth’s atmosphere as a shielding mass, at sea level of about 1033 g of matter per cm3. Indeed high energy protons will generate an atmospheric nuclear cascade, and these high-energy cascading particles have enough energy to continue the process to sea level. The cycles are then more entertained more transparent near the top of the atmosphere reflecting more solar activity and not masked nor dilated with other phenomena. The decrease of T1 with latitude is related to the configuration of the Earth’s magnetic field having a kind of shield less important on higher latitudes which are then more sensitive to solar phenomena particles and activity. T2 increases with altitude up to a certain value of latitude around 1 rd (571) and then T2 remains constant versus altitude, T2 decreases with latitude and, it is generally varying between 11.39 and 39.60 years. The T3 cycle increases with altitude with the same shape of variation, T3 increase with latitude parameter remaining in the interval 53.20, 90.93 years. Unlike T1, T2 and T3 seem to have the opposite behaviour, they generally increase with altitude. Indeed the cycles are more dilated with altitude indicating lower sensitivity to certain Earth’s atmospheric phenomena which may be then most intensive on the Earth’s level. However we can notice that T2 has a mixing behaviour between T1 and T3 combining the reflection of solar and Earth’s phenomena. Besides the cycle T3 is longer with higher latitude excluding the Earth’s magnetic field as effect, this may be so purely related to another Earth’s phenomena decreasing or shielding from the cosmic rays intensities and even surmounting the opposite effect of Earth’s magnetic field. We have proceeded in another step to calculate the first three main periods for Tunisia, Egypt, Algeria, Maroc and Libya (Table 2). We have found a first period T1
Fig. 5. Morlet decomposition coefficients for South Africa cases.
ARTICLE IN PRESS N. Zarrouk, R. Bennaceur / Acta Astronautica 66 (2010) 1008–1016
corresponding to the higher decomposition coefficients according to the corresponding training, T1 is ranging from 12 to 18 years for all the studied countries. The
1015
second cycle T2 is varying from 20 to 24 years. The less frequent period T3 is ranging from 44 to 82 years. Thus we illustrate in this first attempt of prediction, the known
Fig. 6. Morlet reconstructed CR variations curves for Tunis, Alger, Cairo, Rabat, and Tripoli.
ARTICLE IN PRESS 1016
N. Zarrouk, R. Bennaceur / Acta Astronautica 66 (2010) 1008–1016
11–12 cycle, with the single exception of Maroc country for all periods T1, T2, and T3 which seems to be overestimated by neural network system and constitute then a singularity. 3.4. Decomposition and reconstruction of cosmic ray variation for virtual stations Morlet reconstruction test: As we have mentioned above we carried out the training test on the totality of stations except one station on which we apply the results of training and confront them after reconstruction of the signal to measured cosmic rays variations of corresponding station. We have repeated this procedure for all examples of stations. Thus these Morlet test reconstructions are used as test as well for suitability of Morlet reconstruction base as for credibility of neural network prediction. In following curves of Fig. 4 we showed Morlet reconstructed tests for cosmic rays variations compared to corresponding measured CR variations for the period of 1982–2000. In the second stage of training we have then used the totality of stations examples in a main training. Once the NN system has converged with a satisfactory precision with respect to original decomposition coefficients, we have determined then the inputs corresponding to Rabat, Alger, Tunis, Tripoli, Cairo geographic coordinates of virtual neutron monitors stations where we want to derive cosmic ray variation. At the end we have presented these inputs to neural network system already trained, we reach then our purpose when we obtain the unknown decomposition coefficients shown in Fig. 5. We have proceeded to reconstruction in order to have the desired variation cosmic ray curves corresponding then to our virtual neutron monitor stations. We have used a four layers network, an input layer containing three neurons corresponding to geographic coordinates, two hidden layers composed of 16 neurons for wR and an output layer composed of 16 neurons. For processing data a supervised learning was applied with a gradient back propagation algorithm. The Morlet reconstructed curves seem to have the same shape with the original curve given (Fig. 6) by Moscow NM, and also similar shapes for the common 20 years period corresponding to the training period between 1982 and 2002, especially having a minimum corresponding to the lowest relative cosmic rays intensities for the period of 1992–1993 and a maximum indicating a peak for the cosmic rays intensities for the five countries for years 1988–1989, a secondary minimum is also appearing for the period 1984–1986. The main 11 years cycle is appearing between 1984 and 1993, the length of the cycle is shorter in this case and it ranges from 7 to 9 years for the North Africa case.
4. Conclusion We have succeeded through this work to extrapolate numerically the network of neutron monitors stations to other points on the Earth’s surface such us Tunisia, Egypt, Algeria, Maroc, Libya and other points which just must be in the training interval that we have used in our neural network system. Through this work we have found the aspect of very known 11-years cycle: T1, we have revealed the variation type of T2 and especially T3 cycles which seem to be induced by particular Earth’s phenomena that will be studied in next works. Thus in this manner we have derived the Morlet decomposition coefficients corresponding to geographic coordinates of different points. Consequently we have plotted relative cosmic rays intensities curves for a certain period suitable to give a good precision in the training phase 16 or 50 years. The main periods modulating galactic cosmic rays variations detected according to Morlet decomposition coefficients and reconstructions for the five virtual stations are around 9–11 years. Neural network in combination with Morlet wavelets seems to construct a very suitable tool for our case of numerical extrapolation from NM stations network. References [1] Moscow neutron monitor /http://helios.izmiran.troitsk.ru/cosray/ main.htmS. [2] /http://cosmicrays.oulu.fi/llS. [3] /http://neutronmonitor.ta3.sk/realtime.php3llS. [4] /http://www.puk.ac.za/physics/data/S. [5] /http://cosray.unibe.ch/S. [6] /http://cgm.iszf.irk.ru/irkt/main.htmS. [7] /http://cr0.izmiran.rssi.ru/clmx/main.htmS. [8] K.O. Manuel, W.B. Ninham, E.S. Friberg, Super fluidity in the solar interior: implications for solar eruptions and climate, Journal of Fusion Energy 21 (2002) 3–4. [9] H. Mavromichalaki, C. Sarlanis, G. Souvatzoglou, S. Tatsis, A. Belov, E. Eroshenko, V.G. Yanke, A. Pchelkin, Athens neutron monitor and its aspects in cosmic ray variations studies, in: Proceedings of the 27th ICRC 2001, 4099, 2001. [10] H. Demuth, M. Beale, M. Hagan, Neural Network Toolbox 5 User’s Guide. [11] V.I. Kozlov, V.V. Markov, Wavelet image of the fine structure of the 11-year cycle based on studying cosmic ray fluctuations during cycles 20–23, Geomagnetism and Aeronomy 47 (2007) 43–51. [12] M.R. Attolini, et al., Planetary and Space Science 23 (1975) 1603. [13] B.F. Chao, I. Naito, Wavelet analysis provides a new tool for studying Earth’s rotation, EOS 76 (1995) 161–165. [14] N. Zarrouk, R. Bennaceur, A wavelet based analysis of cosmic rays modulation, Acta Astronautica 65 (2009) 262–272. [15] J.M. Torres-Moreno, Apprentissage et ge´ne´ralisation par des re´seaux de neurones: e´tude de nouveaux algorithmes constructifs. The se, Institut National Polytechnique de Grenoble, 1997. [16] J. Pe´rez-Peraza, V. Velasco, S. Kavlakov, Wavelet coherence analysis of Atlantic Hurricanes and cosmic rays, Geofisica International 47 (3) (2008) 231–244. [17] J. Pe´rez-Peraza et al., On the trend of Atlantic Hurricane with cosmic rays, in: Proceedings of the 30th International Cosmic Ray Conference, Mexico 2008.