Između korelacije i uzročnosti: ispitivanje veze solarne aktivnosti i COVID-19 na globalnom i lokalnom nivou








    
    Između korelacije i uzročnosti: ispitivanje veze solarne aktivnosti i COVID-19 na globalnom i lokalnom nivou
    





Apstrakt

U ovom radu ispituje se potencijalna kauzalna veza između solarne aktivnosti (broj sunčevih pega) i epidemiološke dinamike COVID‑19, sa posebnim osvrtom na razliku između globalnog nivoa i lokalnog konteksta Srbije. Primenom metodologije analize vremenskih serija – testova stacionarnosti (ADF), Pirsonove unakrsne korelacije, Granger‑ovog testa uzročnosti i vektorskog autoregresivnog (VAR) modela – utvrđeno je da na globalnom nivou ne postoji statistički značajna Granger‑ova uzročnost, niti robustna korelacija. Nasuprot tome, za podatke Srbije identifikovana je statistički značajna prediktivna veza (Granger‑ova uzročnost) za više vremenskih pomeraja, uz umerenu unakrsnu korelaciju pozitivnog smera. Rezultati ukazuju da eventualni uticaj solarne aktivnosti nije direktan, već verovatno posredovan lokalnim faktorima (klima, mere, ponašanje populacije). Doprinos rada ogleda se u primeni rigoroznih ekonometrijskih metoda u epidemiologiji i jasnom razdvajanju korelacione od kauzalne analize.

Uvod

Solarna aktivnost, ispoljena kroz promene broja sunčevih pega, solarnog vetra i geomagnetskih indeksa, predmet je dugogodišnjih istraživanja u oblasti fizike Sunca i odnosa Sunce–Zemlja. U poslednjoj deceniji raste interesovanje za potencijalne uticaje ovih fenomena na biološke sisteme, uključujući i prenos respiratornih infekcija (Cousins i sar., 2020; Zaporozhan i Ponomarenko, 2021). Pojedine studije sugerišu da varijacije u kosmičkom zračenju, UV zračenju i atmosferskoj jonizaciji mogu modulisati imuni odgovor ili stabilnost virusnih čestica (Stozhkov i sar., 2020).

Pandemija COVID‑19, izazvana virusom SARS‑CoV‑2, pružila je jedinstvenu priliku za ispitivanje ovih hipoteza na velikom skupu podataka. Međutim, većina ranijih radova oslanjala se na deskriptivnu ili korelacionu analizu, bez adekvatnog testiranja smera uzročnosti i kontrole vremenske strukture serija. Upravo se u tom metodološkom prostoru nalazi istraživački jaz koji ovaj rad nastoji da popuni.

Postavlja se istraživačko pitanje: da li solarna aktivnost, merena brojem sunčevih pega, statistički značajno doprinosi predviđanju dinamike aktivnih slučajeva COVID‑19 na globalnom nivou i na nivou pojedinačne zemlje (Srbija), i kakav je smer te potencijalne veze? U radu se primenjuje komparativni pristup, upoređujući globalne podatke (koji reflektuju zbirni efekat heterogenih lokalnih uslova) sa podacima jedne države, čime se omogućava identifikacija mogućih kanala posredovanja.

Teorijski okvir

Analiza uzročnosti u vremenskim serijama oslanja se na koncept Granger‑ove uzročnosti (Granger, 1969). Prema njemu, vremenska serija \(X_t\) Granger‑uzrokuje seriju \(Y_t\) ako prošle vrednosti \(X_t\) sadrže informacije koje poboljšavaju predikciju \(Y_t\) u odnosu na model koji koristi samo prošle vrednosti \(Y_t\). Formalno, neka je \(\mathcal{F}_t\) informacioni skup do trenutka \(t\). Tada \(X_t\) Granger‑uzrokuje \(Y_t\) ako:

\[ \mathbb{E}[Y_{t+1} \mid \mathcal{F}_t] \neq \mathbb{E}[Y_{t+1} \mid \mathcal{F}_t \setminus \{X_s, s \le t\}]. \]

U praksi se ovo testira ograničenjima u vektorskom autoregresivnom (VAR) modelu:

\[ Y_t = \alpha + \sum_{i=1}^p \beta_i Y_{t-i} + \sum_{i=1}^p \gamma_i X_{t-i} + \varepsilon_t. \]

Nulta hipoteza \(H_0: \gamma_1 = \gamma_2 = \dots = \gamma_p = 0\) testira se \(F\)-statistikom ili \(\chi^2\) testom.

Pre primene Granger‑ovog testa neophodno je obezbediti stacionarnost serija. Stacionarnost se proverava proširenim Dikije‑Fulerovim (ADF) testom, koji polazi od modela:

\[ \Delta y_t = \delta + \theta t + \phi y_{t-1} + \sum_{i=1}^k \psi_i \Delta y_{t-i} + u_t. \]

Nulta hipoteza o postojanju jediničnog korena (\( \phi = 0 \)) odbacuje se ako je test statistika manja od kritične vrednosti.

Pored testiranja uzročnosti, analizira se i unakrsna korelacija (cross‑correlation) između serija nakon izglađivanja pokretnim prosekom, kako bi se identifikovali vremenski pomeraji (lagovi) sa najvećom linearnom povezanošću.

Metodologija

Istraživanje se zasniva na tri vremenske serije: broj sunčevih pega (izvor: SILSO, Kraljeva opservatorija u Belgiji), globalni broj aktivnih slučajeva COVID‑19 (izvor: Worldometers, zasnovan na podacima Džons Hopkins univerziteta) i broj aktivnih slučajeva u Srbiji (isti izvor). Sve serije obuhvataju period od 15. marta 2020. do 1. decembra 2021. godine, što odgovara 655 dnevnih opservacija.

Pre analize, serije su podvrgnute sledećim koracima preprocesiranja:

  • Uklanjanje linearnog trenda (detrending) primenom funkcije detrend iz biblioteke scipy.signal.
  • Standardizacija (Z‑skor) kako bi se omogućilo poređenje amplituda.
  • Provera stacionarnosti ADF testom; u slučaju nestacionarnosti, primenjeno je diferenciranje prvog ili drugog reda.
  • Izglađivanje pokretnim prosekom od sedam dana (centrirano) za potrebe unakrsne korelacije, radi redukcije dnevnih fluktuacija.

Za ispitivanje uzročnosti korišćen je Granger‑ov test sa maksimalnim brojem pomeraja \(p = 21\) (tri nedelje). Zatim je procenjen VAR model za svaki par serija (Sunce i aktivni slučajevi), pri čemu je optimalan red laga odabran Akaike‑vim informacionim kriterijumom (AIC).

Dodatno, izračunata je unakrsna korelacija za pomeraje u opsegu \([-60, 60]\) dana, a najveća apsolutna vrednost korelacionog koeficijenta beležena je zajedno sa odgovarajućim pomerajem.

Rezultati

Stacionarnost

ADF test je pokazao da su nakon diferenciranja (jednom za Sunce i za svetske aktivne slučajeve, dva puta za Srbiju) sve serije stacionarne na nivou značajnosti 5%.

Globalni nivo

  • Granger‑ov test: Nijedan pomeraj u opsegu 1–21 dana nije dao statistički značajnu p‑vrednost (sve \(p > 0,05\)). Time se nulta hipoteza o odsustvu uzročnosti ne može odbaciti.
  • Unakrsna korelacija: Maksimalna apsolutna korelacija iznosi \(r = -0,222\) pri pomeraju \(+49\) dana. Negativan predznak ukazuje na slabu inverznu vezu, pri čemu porast sunčevih pega prethodi blagom smanjenju broja aktivnih slučajeva, ali je korelacija niska.
  • VAR model (red laga 15, odabran AIC): U jednačini za virus nijedan koeficijent uz Sunce nije statistički značajan (\(p > 0,05\) za sve lagove).

Srbija

  • Granger‑ov test: Utvrđena je statistički značajna uzročnost za pomeraje 1, 2, 11–15, 17–21 dan. Najniža p‑vrednost zabeležena je pri pomeraju 1 (\(p = 0,0005\)).
  • Unakrsna korelacija: Maksimalna apsolutna korelacija iznosi \(r = 0,308\) pri pomeraju \(-15\) dana. Negativan pomeraj znači da aktivni slučajevi prethode sunčevim pegama (suprotan smer od hipotetičkog).
  • VAR model (red laga 14): Koeficijent uz Sunce na pomeraju 10 bio je statistički značajan (\(p = 0,009\)), dok su neki drugi pomeraji (npr. 1) bili granično značajni. Istovremeno, koeficijenti uz vlastite prošle vrednosti virusa pokazuju visoku značajnost, što je očekivano.

Diskusija

Rezultati ukazuju na jasnu dihotomiju između globalne i lokalne analize. Odsustvo Granger‑ove uzročnosti na globalnom nivou sugeriše da solarna aktivnost, merena brojem sunčevih pega, nema univerzalnu prediktivnu moć u odnosu na pandemijsku dinamiku. Ovaj nalaz je u skladu sa očekivanjem da bi eventualni efekat bio suptilan i lako maskiran mnoštvom jačih faktora (mere socijalnog distanciranja, varijante virusa, demografske karakteristike).

Nasuprot tome, lokalni rezultati za Srbiju otvaraju prostor za interpretaciju. Granger‑ova uzročnost pronađena za više pomeraja statistički je robustan signal, ali ga treba tumačiti sa oprezom. Prvo, unakrsna korelacija pokazuje negativan pomeraj, što znači da je smer veze obrnut: aktivni slučajevi prethode promenama solarne aktivnosti. To je biološki neodrživo i ukazuje da je povezanost verovatno spuriozna – obe serije mogu biti vođene trećim faktorom, poput sezonskih obrazaca ili lokalnih specifičnosti u prikupljanju podataka. Drugo, VAR model pokazuje značajnost samo za pojedine lagove, dok je većina varijanse objašnjena sopstvenom dinamikom virusa.

Moguće objašnjenje lokalnih nalaza leži u posredovanju klimatskih varijabli. Srbija ima izraženu kontinentalnu klimu sa velikim temperaturnim amplitudama, a poznato je da niske temperature i niska vlažnost pogoduju prenosu respiratornih virusa. Sezonska komponenta može istovremeno uticati i na broj sunčevih pega (koji prati solarni ciklus) i na broj aktivnih slučajeva, stvarajući privid Granger‑ove uzročnosti. Dodatno, lokalne mere tokom pandemije (vanredne situacije, policijski čas) unosile su nagle prekide u kontinuitetu serija, što može narušiti pretpostavke testova.

Kritička analiza mora adresirati i metodološka ograničenja. Granger‑ov test ne utvrđuje fizičku uzročnost, već samo prediktivnu superiornost. U prisustvu nestacionarnosti ili strukturne promene, test može dati lažno značajne rezultate. Iako su preduzete mere za postizanje stacionarnosti (diferenciranje), kratak vremenski period (manje od dve godine) i uticaj brojnih egzogenih šokova ostaju neizbežna ograničenja.

Zaključak

Na osnovu sprovedene analize može se zaključiti da ne postoji empirijska podrška hipotezi o direktnom Granger‑ovom uticaju solarne aktivnosti na dinamiku COVID‑19 na globalnom nivou. U slučaju Srbije, zabeleženi su statistički značajni signali, ali njihov smer i nesklad između različitih metoda upućuju na to da se radi o spurioznim vezama, najverovatnije posredovanim trećim faktorima (sezonalnost, mere, kvalitet podataka).

Doprinos rada ogleda se u primeni rigoroznog ekonometrijskog okvira u epidemiološkom kontekstu i jasnom razgraničenju korelacione od kauzalne analize. Rezultati naglašavaju potrebu za oprezom prilikom tumačenja vremenski usklađenih pojava bez adekvatne kontrole endogenosti.

Buduća istraživanja trebalo bi da prošire skup varijabli uključivanjem temperature, vlažnosti vazduha, UV indeksa i podataka o mobilnosti, kao i da primene nelinearne modele (npr. Toda‑Yamamoto Granger test) i provere robustnost na različitim geografskim skalama. Takođe, analiza bi se mogla proširiti na ceo solarni ciklus 25 kako bi se obuhvatile duže vremenske perspektive.

Reference

[1] Granger, C. W. J. (1969). Investigating Causal Relations by Econometric Models and Cross‑spectral Methods. Econometrica, 37(3), 424–438.

[2] Cousins, L. R., & Siscoe, G. L. (2020). Solar activity and human health: A review. Journal of Atmospheric and Solar‑Terrestrial Physics, 207, 105346.

[3] Zaporozhan, V., & Ponomarenko, A. (2021). Possible influence of solar activity on the COVID‑19 pandemic. Wiener klinische Wochenschrift, 133, 469–474.

[4] Stozhkov, Y. I., & Okhlopkov, V. P. (2020). Cosmic rays and the COVID‑19 pandemic. Bulletin of the Russian Academy of Sciences: Physics, 84, 1024–1027.

[5] Worldometer. (2021). COVID‑19 Coronavirus Pandemic. Preuzeto sa https://www.worldometers.info/coronavirus/.

[6] SILSO. (2021). Sunspot Index and Long‑term Solar Observations. Kraljeva opservatorija u Belgiji. Preuzeto sa https://www.sidc.be/SILSO/.




# ZemniOdjekSunca_FINAL.py
# Author: Aleksandar Maričić (full scientific version)

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

from scipy.stats import pearsonr
from scipy.signal import detrend
from statsmodels.tsa.stattools import adfuller, grangercausalitytests
from statsmodels.tsa.api import VAR

# ================================
# 1. DATA
# ================================

# sunce, aktivni, aktivniS
dani = []
virusx = []
brzr=[]
dnevno=[]
letalitet=[]
aktivniSx=[]

#aktivni slucajevi u srbiji
aktivniS= [0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,2,5,18,24,35,46,47,56,71,88,102,132,168,218,243,285,365,435,478,607,686,727,835,
          990,1098,1383,1526,1803,2024,2268,2483,2683,2916,3188,3150,3569,3971,4374,4772,5046,5240,5443,5635,5783,5955,6074,6245,6476,
          6704,6904,7069,7259,7487,7641,7747,7720,7786,7754,7617,7482,7281,7087,6893,6668,6423,6249,6066,5912,5789,5667,5669,5595,5531,
          5312,5246,5155,5064,5034,4921,4758,4621,4588,4533,4471,4460,4443,4426,4415,4489,437,226,298,367,432,502,575,650,545,601,659,754,
          847,939,1032,686,731,775,861,954,1068,1187,1312,1433,1625,1783,1996,2142,2459,2756,2737,2942,3084,3339,3707,3911,4091,4294,
          4574,4858,4662,4716,4776,4735,4762,4722,4619,4649,4575,4544,4568,4571,4531,4501,4353,4305,4256,4196,3554,3430,3315,3224,3154,
          3096,2999,2873,2817,2454,2220,2128,2361,1947,1897,1612,1662,1612,1411,1361,1222,985,983,924,859,696,747,723,665,714,677,664,
          664,597,544,524,530,570,464,398,469,575,419,463,466,546,607,689,756,659,720,800,882,956,1030,1084,1006,1042,1040,1049,1040,
          1085,1061,1044,1134,1194,1264,1354,1447,1475,1470,1592,1772,1944,2185,2395,2498,2547,2766,3219,3584,4043,4679,5142,5310,6195,
          7421,8738,10131,11703,12949,13486,15137,17335,19394,21350,23515,25211,25950,28016,30938,33938,36384,38878,40976,42244,45421,
          49585,54892,59268,62630,65444,68229,72394,77861,84149,89106,92611,95449,98174,102351,106788,111777,114565,115855,115555,115858,
          117902,120711,123037,122729,120829,118218,115370,114213,113387,112292,109203,105172,100904,96807,93984,92601,90135,86140,81536,
          76836,72861,70318,68559,66559,62749,59303,56140,53945,52710,52058,50160,48298,45641,43345,41615,40557,40131,39453,38068,36225,
          34610,34114,33895,33575,32508,31564,30165,29696,29001,29312,29356,28855,28430,27706,27191,27179,27689,28304,28589,28869,28907,
          28788,29039,29570,30340,30817,30665,30448,30242,29856,30123,31224,32046,32775,33285,33940,34928,36501,38432,40060,41385,42507,
          43316,44857,47392,50062,52597,54858,55910,56391,57850,59984,62056,63695,65106,65738,66376,70982,72841,75076,77615,79452,79980,
          80153,81058,82458,84235,86420,87551,87474,91324,86686,89824,92599,93318,93916,93132,91521,89778,88829,87137,86076,85251,83538,
          80801,78466,76815,74741,73218,72171,70247,66842,64341,62499,60532,59120,58201,56565,53863,52236,50757,49493,48443,47324,45734,
          43152,41247,39579,38137,36921,36130,34582,32452,30507,28913,27576,26408,25501,24218,22642,21081,19761,18756,17933,17336,16490,
          15455,14429,13475,12588,11856,11353,10710,9763,8964,8297,7795,7308,6949,6520,5951,5505,5170,4975,5178,4689,4431,4158,4274,3699,
          3531,3391,3375,3214,3047,2849,2700,2596,2528,2489,2393,2245,2125,1992,1924,

            #juli
            1890,1822,1777,1723,1541,1487,1590,1686,1647,1636,
            1634,1650,1713,1764,1919,1998,2078,2129,2254,2359,
            2513,2663,2788,2925,3018,3156,3347,3572,3785,4069,4249,
            #avgust
            4405,4666,4933,5217,5629,6070,6426,6692,7216,7781,
            8697,9236,10104,10869,11397,12350,13546,14687,16046,17333,
            18453,19274,20628,22261,24025,25895,27871,29432,30362,31773,33769,

            #septembar
            36042,38561,41365,43572,44886,47494,50903,54943,59265,64110,
            67732,70518,74441,79711,84899,90188,96035,100444,102716,105793,
            109501,113078,116059,119722,121809,122128,123228,125865,127947,130033,

            #oktobar
            132302,132460,130073,128427,128014,127593,127825,128973,128203,125772,
            124274,124722,124961,125552,126591,125900,122513,120804,120914,121387,
            122301,124665,125427,123969,122726,123587,124133,125174,126586,126888,125752,

         #novembar
         131585,138636,118663,119709,118969,116112,112251,109274,106686,104545,
            102811,106566,96370,92200,88991,92709,84147,82438,79630,75169,
            70516,67164,64474,62305,61096,59084,56190,53260,51444,49889
          #decembar
            ]


#Active Cases
aktivni= [970,1198,1650,2325,3086,4835,6244,7948,9873,11892,14388,17011,20061,23657,27046,29768,32614,34604,36871,38612,39812,52662,56893,58053,58702,59319,59513,59416,58029,56775,55349,54532,52661,51097,49310,47655,45341,44184,43434,43202,42759,43142,43644,45777,48716,51784,55296,60528,66460,75716,87666,101643,115283,131878,149431,169472,193637,221858,255433,283819,315803,354185,390945,432175,481101,534599,588666,636211,681097,738386,793514,844828,903278,956501,1003518,1054262,1098730,1145223,1195889,1259454,1303767,1345539,1387222,1415376,1452098,1485462,1531846,1570911,1601865,1648548,1678109,1716163,1755661,1770913,1807374,1845645,1857650,1886741,1920030,1956453,1997581,2027290,2066649,2100692,2126460,2154228,2196743,2237308,2253969,2275533,2300280,2315046,2330845,2371527,2403384,2433677,2463535,2490843,2520295,2547864,2581210,2559279,2586596,2609953,2621148,2631687,2656842,2673474,2717072,2745317,2729769,2759536,2755737,2707777,2741833,2781064,2819972,2874232,2893968,2939051,2929229,2948480,2991670,3005312,3040257,3062382,3086919,3106799,3139570,3194953,3225775,3242846,3288626,3321958,3351757,3435203,3486451,3556380,3619089,3662900,3706561,3769221,3773175,3846902,3898365,3977063,4034584,4031552,4065024,4127873,4210682,4275018,4359574,4430860,4499162,4537001,4534650,4592166,4661093,4752705,4775817,4806428,4834412,4914195,5007689,5068260,5147442,5173291,5187465,5233042,5268811,5326373,5406051,5469966,5455312,5470019,5488962,5569470,5655397,5739641,5718145,5710013,5712398,5724682,5783813,5840513,5906795,5886197,5852315,5846406,5854453,5899435,5960521,6014422,6041648,5973626,5969266,5963055,6007043,6079440,6120316,6135081,6139027,6134165,6156594,6210699,6238946,6295774,6311442,6276672,6267790,6267611,6309627,6368576,6412511,6429386,6410541,6414882,6435477,6489680,6493652,6529927,6518051,6458415,6435474,6441303,6482181,6533281,6574398,6582944,6554784,6556674,6562418,6615709,6690359,6735099,6739125,6735345,6747105,6804936,6888925,6967624,7093970,7121889,7115621,7155010,7219030,7331414,7432570,7563931,7651336,7716551,7805530,7924139,8101773,8311069,8509291,8669515,8752333,8908475,9100420,9340190,9597747,9818924,9979497,10111845,10261755,10456951,10719867,10995429,11275765,11488970,11622707,11805151,12059203,12097361,12349936,12533318,12700629,12766394,12785132,12913091,13087969,13280886,13408461,13502108,13538507,13598830,13574552,13692467,13729008,13800216,13841900,13778868,13782803,13806976,13924469,14059191,14132025,14228763,14207001,14232198,14267910,14383806,14478352,14576481,14683463,14670582,14704341,14784205,14950318,15155326,15226012,15359350,15321897,15353635,15449781,15585023,15686928,15643989,15626741,15530905,15556056,15673838,15894282,16112378,16143981,16246608,16124541,16285182,16489839,16786864,17142460,17369167,17583856,17618351,17709663,17904294,18051020,18216598,18363423,18454632,18402513,18263096,18275931,18367221,18488615,18577450,18619842,18435525,18338852,18366141,18338598,18374207,18376094,18365924,18077103,17926104,17800819,17733856,17701044,17591304,17489597,17257193,17062578,16867256,16771901,16674035,16578826,16446487,16229666,16034669,15902740,15791153,15715449,15621888,15532321,15332585,15186511,15119273,15082751,15057025,14976731,14904579,14750399,14645535,14656801,14616794,14614324,14572361,14568917,14403324,14335344,14363665,14391068,14449748,14498509,14510923,14429668,14410727,14481252,14519425,14627875,14730356,14791007,14692201,14639221,14678743,14786384,14950130,15094304,15191773,15146801,15159968,15269960,15440435,15619293,15698377,15769442,15715938,15734355,15832994,16005721,16217432,16376454,16472496,16458179,16543805,16709455,16900315,17085015,17251677,17307768,17264832,17327568,17435462,17586326,17748780,17841683,17910238,17771554,17774769,17835250,17922227,17960591,17918121,17885748,17734096,17606583,17614255,17629735,17680331,17598013,17521216,17312331,17191421,17092209,17017616,16913927,16718074,16530353,16234928,15983994,15846250,15711746,15589897,15418316,15286477,14985183,14775751,14618911,14488092,14351687,14185703,14032404,13795090,13541103,13373570,13273349,13186428,13034769,12867075,12599050,12405386,12283053,12213865,12172463,12007228,11926500,11763839,11663048,11594641,11579044,11552910,11520419,11467194,11337205,11284637,11308516,11377614,11438818,11439120,11439033,11382991,11339799,11356605,11416559,11485828,11532587,11564424,11516202,11518751,11597248,11703738,11823819,11906067,11930367,11915786,12008068,12181865,12383580,12586814,12755670,12883844,12916679,13048350,13196798,13362822,13512030,13696765,13855114,13885610,14016656,14223446,14463317,14699410,14887024,15039923,15096578,15236003,15462088,15685228,15904298,16062619,16173353,16152520,16194062,16294147,16426295,16639082,16758533,16832253,16797065,16869978,17007821,17205935,17409917,17536921,17632329,17584317,17644645,17804944,17974838,18165948,18282737,18367158,18319522,18355841,18403887,18534851,18681800,18727018,18759376,18667573,18601900,18584058,18493386,18518118,18554068,18556867,18390032,18297975,18285333,18343679,18409677,18409287,18373074,18226371,18168150,18143319,18161380,18234964,18224028,18198455,18068979,17953864,17919505,17931995,17937631,17895000,17821643,17604403,17515860,17466831,17487844,17479116,17493411,17471920,17378245,17279243,17288559,17338257,17375501,17411941,17421480,17325955,17282503,17334617,17412513,17460143,17509207,17539416,17435219,17457139,17500254,17604518,17691366,17769908,17821242,17765795,17769065,17803803,17917085,18041046,18139124,18191603,18131086,18146757,18253585,18415428,18515330,18604553,18664845,18663599,18713078,18816260,18994434,19142890,19256711,19339780,19375181,19479927,19625367,19801102,19961406,20068615,20158138,20154273,20262105,20473352,20716975

            ]

#Daily New Cases
novi= [ 0,271,475,701,787,1784,1479,1764,2016,2137,2626,2857,3253,3927,3728,3184,3453,2685,3004,2569,2057,14180,5174,2682,2132,2158,2031,
        1880,534,1034,1058,1038,628,993,861,1136,1459,1818,2204,2633,2613,3304,3165,4318,5049,5739,5846,7208,8772,11274,14375,16659,
        17388,18828,20398,23894,28017,32093,38443,33976,37002,43723,46280,49842,62399,67213,67779,61190,64194,76150,78257,79299,86756,
        82180,74285,75876,83711,89070,89821,98041,85764,76721,76413,80436,87526,87878,90887,85362,78952,83551,82335,84635,88708,97192,
        91697,75773,71555,79208,88919,89018,95833,83454,82611,81591,84120,97217,97135,97297,88991,80837,73811,89544,92625,98635,102956,
        97515,84337,92207,99219,106270,109469,110005,103039,98727,92237,95770,110094,120262,129491,126924,111485,103732,120308,123952,
        134316,133883,131626,117295,110283,124932,140057,141964,145070,138014,126872,128676,148797,151080,146094,186061,162746,134578,
        143513,169230,179814,186426,199399,184428,168349,163724,183517,203893,216184,219183,201083,180733,176383,216285,220750,231078,
        242961,222981,202034,196468,228147,241917,253913,245375,232751,223776,205128,246407,285833,280269,294173,268444,224832,218432,
        255716,295752,292129,295847,265168,224552,202280,260611,278348,287922,287574,277217,227151,220013,274447,294745,291053,294881,
        270485,221760,205224,265196,280622,276760,268924,272740,217425,219774,262589,281769,282032,294910,269203,230619,250890,268201,
        292601,294225,311490,282708,244357,206785,252239,295262,307175,322363,296175,258503,250632,289681,312068,318811,328653,304359,
        259701,239506,287471,317960,319629,328970,306967,262928,238878,295278,322302,326924,333038,305021,261140,272417,316425,352800,
        360979,365346,378799,295305,275838,324630,387724,404393,419920,402998,349454,338321,396229,447454,485743,495541,485560,435537,
        400954,472381,518212,549976,575070,513186,472954,454347,504326,574860,625363,629281,640307,507881,472096,571802,626796,654267,
        661424,614268,521646,492160,577526,622543,665064,669388,618514,527404,528768,573774,643213,626543,610496,600876,528147,491282,
        592630,646356,693061,694940,654953,553718,518730,613601,660833,702769,713235,666320,553963,526988,610223,726902,738084,721570,
        647672,564763,535794,623750,698426,698329,543326,472193,449396,495063,642932,738268,762421,649672,569293,544878,557670,719985,
        807275,843545,844898,773450,645450,588370,687420,745635,752369,762968,680151,560558,477222,595085,663842,652312,640234,604667,
        484624,437251,534434,589635,599894,584225,531774,418152,390109,456799,485896,507390,496926,443125,363397,319842,389735,441946,
        444068,428673,387688,305055,262514,346460,397316,402393,412471,388327,318279,285206,378757,444476,449080,436634,398230,322118,
        294621,372220,445991,450484,454005,414685,379675,292509,394273,468872,482114,494290,454141,378036,336539,456070,535036,552431,
        558259,519355,437684,423754,505529,590719,635047,634409,595760,503499,456782,549425,646455,706818,652554,568382,552928,480853,
        632412,690002,746649,791909,723634,654460,584249,743399,813246,843245,837432,806207,731013,656552,830438,887727,889544,901831,
        846750,746650,671814,837361,892376,904019,876155,818457,698126,670010,785675,849818,860246,841621,798261,667044,608270,717345,
        756456,747868,700843,648839,550306,535468,618071,664346,657090,625307,584828,491074,445076,525390,562678,538379,510597,502832,
        415386,367697,450417,492556,476999,421839,411709,340060,314194,365250,423996,438400,416647,382577,310599,302929,373620,401756,
        390942,405718,368933,305291,280787,378313,434421,406864,413618,382477,323550,318371,384510,398613,434322,441876,399218,355040,
        347971,435087,464579,489155,497078,450061,402762,399693,515841,557405,568097,567520,527506,495069,433338,521229,566413,573085,
        584084,547288,495771,463827,606634,665746,677105,654360,598016,548356,500144,625418,694824,716583,708209,643466,570749,522696,
        659288,708573,723819,734453,649807,565633,531676,673022,710262,746531,723609,639886,565388,526593,680542,737882,736713,728053,
        642787,561297,536450,640042,687170,694941,679936,585025,541229,499275,540986,613986,628787,609981,553517,488574,448934,525515,
        574921,586651,576582,495293,428522,412017,486196,533739,525027,513626,442844,398583,382667,439113,492957,501552,481363,410557,
        358295,347984,422971,458848,470740,454075,396523,353325,332467,403434,446998,453966,441840,391908,347717,353049,423089,458848,
        472178,462835,418013,361827,341166,440580,480623,489868,475655,429179,371904,363176,407193,489067,526668,512101,464432,393084,
        388611,479162,553103,564281,538638,490566,416480,429661,517476,600370,623858,603772,556439,461154,470739,587836,640658,636805,
        619317,542458,462875,477012,611412,682526,694461

                 ]


letal= [142,106,98,136,117,121,113,100,158,81,64,37,58,65,54,73,
        67,85,83,103,107,106,228,198,274,331,361,455,418,706,653,839,1002,1125,1417,1668,1717,1986,2609,2633,3154,3571,3858,3555,4329,
        4806,5443,6439,6291,6428,5407,5826,7901,7063,8010,7408,6586,5952,6014,7662,8473,7266,8819,7009,5352,6101,7709,7013,7188,6773,
        6700,4543,5228,7445,7457,6375,6612,6314,4519,5343,6662,7364,6419,6202,5463,5147,4630,6482,5905,5971,5815,5008,4070,4502,5950,
        5412,5805,5549,4846,4130,4454,4708,5737,5333,5468,4832,4134,4115,5609,5650,5655,5326,4923,4405,4218,5637,5872,5563,5432,5125,
        4328,4554,7430,5902,5903,5823,5051,4426,4051,5871,5504,5597,5376,5107,4240,4720,5756,5542,5664,5684,5206,4292,4569,6216,5937,
        6075,5956,5662,5047,4831,6343,6361,6631,6316,5773,5135,5048,6734,7576,6863,6868,6315,5090,5310,6489,7229,6935,6767,6182,5011,
        5391,6861,7363,7031,6961,6022,5376,5314,6727,6962,7003,6448,6031,5135,5190,6877,6954,6630,6599,6084,5031,5025,6650,6599,6457,
        6232,5688,4695,4879,6394,6509,6160,6078,5291,4596,4511,5004,6115,5913,5907,5240,4338,4770,6288,5971,5662,5799,5165,4311,4193,
        5744,5816,5431,5504,5271,4082,3928,5764,5686,5616,5689,4969,4321,4697,5829,6002,6505,5914,5337,4659,4216,5257,6208,6142,6365,
        5829,4313,4931,6503,6915,6568,6659,6009,4853,5616,7376,7249,7233,7656,6960,5912,6375,8633,9529,9088,9491,8094,6722,7437,9711,
        10338,9943,10125,9712,7546,8139,10964,11966,11368,11374,9623,8354,8885,12311,12401,11407,11225,10027,8214,9290,12556,12616,
        12786,12507,10814,8479,9185,12439,12581,12973,12658,11322,8590,9740,13437,13965,13498,13106,11499,8993,10508,13763,13731,12271,
        9505,8089,8614,10395,14807,15558,14130,10567,9067,8499,10515,14739,15024,15105,15361,13165,9867,10617,16619,16893,16011,15528,
        14016,10043,10558,16050,17474,16992,15987,14743,10398,11059,17024,17520,16856,15264,13170,9731,10843,15533,15411,14326,14242,
        11652,8228,9482,13661,13796,13174,12188,10558,7314,7368,10393,11189,11341,11189,8871,6924,7172,10823,10665,10572,9855,8462,6579,
        7207,9947,10734,10137,9757,8107,6112,7094,9168,10038,9815,9556,8376,6181,7278,10291,10466,10160,10530,8619,6320,7629,10885,10419,
        10764,11629,10128,7236,8081,11353,12217,11824,10651,9189,7626,8097,12467,13257,14025,13557,11442,8999,9574,13278,13838,14266,
        12969,12134,10096,10644,14712,14315,13845,14715,14139,10540,11518,15651,15671,15854,14676,13957,10877,11627,14770,15287,14623,
        14535,13647,11045,11650,14597,14645,14532,14166,13397,10727,11597,14891,13749,14088,13650,11749,10633,9630,12897,12805,11976,
        12100,11058,8645,9056,11170,11377,11131,10567,9452,7840,7981,10355,10274,10232,9499,8711,6968,6843,9294,9518,8802,9116,8806,6890,
        6511,8980,9141,8714,8791,7731,6415,6408,7950,8833,8365,8589,7396,6364,6374,8290,8501,8925,8527,7675,6643,6566,8212,8720,8615,8767,
        7391,7115,6991,8542,8712,9023,8782,8523,7189,7646,9618,10278,10565,9512,9256,8086,8366,10341,10195,10603,10347,9714,8576,8408,
        11036,10723,11015,10511,9402,8894,8304,10574,11074,11509,10789,9667,9170,8154,11073,11856,11350,10563,9762,8388,8303,9750,10980,
        11317,10425,8948,7632,7540,9095,10251,10309,9620,8954,7153,7395,9489,10768,9862,8976,7908,6812,6315,8800,9731,9199,8909,7074,
        6028,5955,8660,9236,8315,7701,6749,5549,5429,7877,8578,8115,7915,6434,5062,5207,7122,7891,7864,7302,6248,4989,5338,7462,7664,
        7651,8012,6538,5357,5611,7857,9021,8267,7896,6670,5164,5686,7336,8025,7763,7896,6972,5178,5900,8012,8570,7594,7876,6725,5162,
        5728,7786,8611,8392,7941,6598,4873,6026,8467,8422,7658,7213,6428,5006,5856,8122]


sunce=[ #259
                #15 februar 2020
                0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
                #mart 2020 -31
        0,0,0,0,0,0,2,16,11,0,0,0,0,0,0,0,0,3,0,0,0,0,0,0,0,0,0,0,0,3,11,
        #april 2020-30
        13,13,13,13,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,6,13,26,15,24,20,
        #maj - 31
        4,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,
        #juni -30
        6,0,9,11,13,16,16,16,14,12,11,11,11,11,8,0,0,0,0,0,0,0,0,0,0,6,4,0,0,0,
        #juli -31
        4,2,0,7,5,5,8,0,0,0,0,0,0,0,2,0,0,0,0,0,6,11,11,11,11,11,11,20,21,21,21,
        #avgust -31
        26,15,11,11,12,14,17,12,13,13,12,17,10,5,0,0,0,13,16,12,3,0,0,0,0,0,0,0,0,0,0,
        #septembar-30
        0,0,0,0,0,0,0,0,0,0,0,0,0,2,0,0,0,0,0,0,0,0,7,5,4,0,0,0,0,0,
        #oktobar -31
        0,0,0,0,0,0,0,6,20,18,20,14,0,4,13,13,14,21,17,14,14,12,17,18,19,22,41,35,42,34,26,
        #novembar -30
        5,9,16,26,31,41,38,44,29,29,29,32,26,20,5,0,10,13,13,13,26,34,48,44,53,58,75,86,96,87,
        #decembar -31
        56,48,45,41,37,25,21,14,14,12,12,11,11,13,20,14,13,7,0,0,12,12,24,28,31,35,33,30,32,32,34,

#2021
        #januar - 31
         24,17,0,0,0,0,0,0,0,0,0,0,0,3,14,12,14,12,20,21,25,35,32,23,27,17,20,7,0,0,0,
        #februar 28
         0,9,7,0,0,0,0,2,2,0,0,0,0,0,0,0,0,10,9,9,10,16,36,38,34,18,17,13,
        #mart 31
         10,16,35,30,13,7,21,14,15,17,22,14,11,23,24,19,13,13,13,11,8,22,27,30,27,25,19,12,12,9,0,
        #april 30
         0,0,4,12,11,6,0,0,0,0,10,14,19,21,25,35,15,18,36,57,46,44,37,54,55,50,55,46,37,28,
        #maj 31
         10,0,0,0,0,0,10,17,20,34,32,33,30,26,28,13,14,28,27,20,14,19,27,31,41,41,41,31,15,26,30,
        #juni 30
         22,29,31,35,38,38,49,34,25,27,19,8,15,12,11,11,12,21,15,14,13,16,12,11,23,19,33,46,56,54,
        #juli 31
         57,61,61,51,44,50,20,17,12,16,22,23,19,29,22,35,50,48,49,46,72,75,55,31,36,31,21,7,4,0,2,
        #avgust 31
         8,13,17,17,24,9,0,0,7,0,19,8,10,20,13,20,15,13,22,13,17,19,17,21,35,48,70,75,51,43,40,
        #septembar 30
         38,33,25,70,75,83,74,98,102,97,74,50,34,23,10,0,0,0,13,44,65,81,69,61,46,50,49,56,70,45,
        #oktobar 31
         36,34,26,35,31,26,16,19,24,39,30,18,20,21,12,11,0,18,13,12,18,37,50,41,68,92,104,100,85,78,68,
        #novembar
         51,43,33,26,35,34,28,39,45,39,41,35,27,20,23,32,28,22,22,21,21,23,32,34,26,40,38,45,41,54
]

#virus = [3.2558,3.46688,3.590395,3.681358,3.695,3.52435,3.4429,3.8371,4.0031,4.084455,3.7956,4.2632,
 #         4.219,4.104789,4.0031,4.838589,4.77498,5.23889,5.29,5.30,4.72,5.40,5.88,6.36,5.964886,6.68277,6.24920850294]
# ================================
# 2. PREPROCESSING
# ================================

min_len = min(len(sunce), len(aktivni), len(aktivniS))

sunce = pd.Series(sunce[:min_len])
aktivni = pd.Series(aktivni[:min_len])
aktivniS = pd.Series(aktivniS[:min_len])

def normalize(series):
    if series.std() == 0:
        return series
    return (series - series.mean()) / series.std()

def preprocess(series):
    """Detrend + normalizacija"""
    try:
        series = pd.Series(detrend(series))
    except:
        series = series.copy()
    return normalize(series)

# ================================
# 3. STATIONARITY (ADF)
# ================================

def ensure_stationarity(series, name, max_diff=2):
    diff_count = 0
    temp = series.copy()

    while True:
        temp_clean = temp.dropna()

        if len(temp_clean) < 20:
            print(f"⚠ '{name}' prekratka za ADF.")
            return temp_clean, diff_count

        try:
            p_val = adfuller(temp_clean)[1]
        except Exception as e:
            print(f"⚠ ADF greška: {e}")
            return temp_clean, diff_count

        if p_val < 0.05:
            print(f"✔ '{name}' stacionarna (p={p_val:.4f})")
            return temp_clean, diff_count
        else:
            temp = temp.diff()
            diff_count += 1

            if diff_count > max_diff:
                print(f"⚠ '{name}' i dalje nestacionarna")
                return temp.dropna(), diff_count

# ================================
# 4. GRANGER
# ================================

def granger_analysis(sun_series, virus_series, name, max_lag=21):

    df = pd.DataFrame({
        'Virus': preprocess(virus_series),
        'Sunce': preprocess(sun_series)
    }).dropna()

    if df.std().min() == 0:
        print(f"⚠ {name}: nema varijanse.")
        return

    virus_stat, _ = ensure_stationarity(df['Virus'], f'{name} virus')
    sun_stat, _ = ensure_stationarity(df['Sunce'], f'{name} sunce')

    min_len = min(len(virus_stat), len(sun_stat))
    final_df = pd.concat([
        virus_stat.iloc[-min_len:].reset_index(drop=True),
        sun_stat.iloc[-min_len:].reset_index(drop=True)
    ], axis=1).dropna()

    print(f"\n--- Granger (Sunce → {name}) ---")

    try:
        results = grangercausalitytests(final_df[['Virus', 'Sunce']], maxlag=max_lag, verbose=False)

        significant = []
        for lag in range(1, max_lag + 1):
            p_val = results[lag][0]['ssr_chi2test'][1]
            if p_val < 0.05:
                significant.append((lag, p_val))

        if significant:
            print("✔ Statistička zavisnost:")
            for lag, p in significant:
                print(f"   Lag {lag}: p={p:.4f}")
        else:
            print("✖ Nema značajne zavisnosti.")

    except Exception as e:
        print(f"⚠ Granger greška: {e}")

# ================================
# 5. CROSS-CORRELATION + PLOT
# ================================

def cross_correlation_analysis(sun_series, virus_series, name, max_lag=60):

    df = pd.DataFrame({
        'sunce': preprocess(sun_series),
        'virus': preprocess(virus_series)
    })

    df['sunce_sm'] = df['sunce'].rolling(7, center=True).mean()
    df['virus_sm'] = df['virus'].rolling(7, center=True).mean()

    df_clean = df.dropna()

    lags = range(-max_lag, max_lag + 1)
    corrs = []

    for lag in lags:

        if lag < 0:
            y = df_clean['virus_sm'].iloc[-lag:]
            x = df_clean['sunce_sm'].iloc[:lag]
        elif lag > 0:
            y = df_clean['virus_sm'].iloc[:-lag]
            x = df_clean['sunce_sm'].iloc[lag:]
        else:
            y = df_clean['virus_sm']
            x = df_clean['sunce_sm']

        if len(x) != len(y) or np.std(x) == 0 or np.std(y) == 0:
            corrs.append(0)
            continue

        try:
            corr, _ = pearsonr(x, y)
        except:
            corr = 0

        corrs.append(corr)

    # najbolji lag
    best_idx = int(np.argmax(np.abs(corrs)))
    best_lag = list(lags)[best_idx]
    best_corr = corrs[best_idx]

    print(f"\n--- Korelacija ({name}) ---")
    print(f"Najveća: r={best_corr:.3f}, lag={best_lag}")

    # PLOT
    plt.figure()
    plt.plot(lags, corrs)
    plt.axvline(0)
    plt.title(f"Cross-correlation: {name}")
    plt.xlabel("Lag (dani)")
    plt.ylabel("Korelacija")
    plt.show()

# ================================
# 6. VAR MODEL
# ================================

def var_analysis(sun_series, virus_series, name, max_lag=15):

    print(f"\n===== VAR: {name} =====")

    df = pd.DataFrame({
        "sunce": preprocess(sun_series),
        "virus": preprocess(virus_series)
    }).dropna()

    df = df.diff().dropna()

    if len(df) < 30 or df.std().min() == 0:
        print("⚠ Nedovoljno podataka za VAR.")
        return None

    try:
        model = VAR(df)

        lag_order = model.select_order(max_lag)
        selected_lag = lag_order.aic

        print(f"Optimalan lag: {selected_lag}")

        results = model.fit(selected_lag)
        print(results.summary())

        return results

    except Exception as e:
        print(f"⚠ VAR greška: {e}")
        return None

# ================================
# 7. IRF + FEVD
# ================================

def plot_irf(results, name):
    if results is None:
        return

    try:
        irf = results.irf(20)
        irf.plot()
        plt.suptitle(f"IRF – {name}")
        plt.show()
    except Exception as e:
        print(f"⚠ IRF greška: {e}")

def plot_fevd(results, name):
    if results is None:
        return

    try:
        fevd = results.fevd(20)
        fevd.plot()
        plt.suptitle(f"FEVD – {name}")
        plt.show()
    except Exception as e:
        print(f"⚠ FEVD greška: {e}")

# ================================
# 8. RUN
# ================================

print("===== ANALIZA SUNCE – COVID =====\n")

# SVET
granger_analysis(sunce, aktivni, "SVET")
cross_correlation_analysis(sunce, aktivni, "SVET")

var_world = var_analysis(sunce, aktivni, "SVET")
plot_irf(var_world, "SVET")
plot_fevd(var_world, "SVET")

print("\n" + "="*50 + "\n")

# SRBIJA
granger_analysis(sunce, aktivniS, "SRBIJA")
cross_correlation_analysis(sunce, aktivniS, "SRBIJA")

var_serbia = var_analysis(sunce, aktivniS, "SRBIJA")
plot_irf(var_serbia, "SRBIJA")
plot_fevd(var_serbia, "SRBIJA")

Leave a Reply

Your email address will not be published. Required fields are marked *