R Under development (unstable) (2024-11-12 r87322 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2024 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> 
> #""""""""""""""""""""""""Detecting periodicity
> #install.packages("expm", repos = "https://cran.r-project.org")
> #install.packages("sn", repos = "https://cran.r-project.org")
> #install.packages("readxl", repos = "https://cran.r-project.org")
> 
> #install.packages("expm")
> #install.packages("sn")
> #install.packages("readxl")
> 
> 
> #library(PerRegMod)
> requireNamespace("PerRegMod")
Loading required namespace: PerRegMod
> requireNamespace("expm")
> requireNamespace("readxl")
> requireNamespace("sn")
> 
> library(PerRegMod)
> library(expm)
Loading required package: Matrix

Attaching package: 'expm'

The following object is masked from 'package:Matrix':

    expm

> library(readxl)
> library(sn)# for skew normal density
Loading required package: stats4

Attaching package: 'sn'

The following object is masked from 'package:stats':

    sd

> 
> 
> 
> n=200
> s=2
> x1=rnorm(n,0,1)
> x2=rnorm(n,0,2)
> x3=rnorm(n,0,3)
> x4=rnorm(n,0,2.7)
> y=rnorm(n,0,2.5)
> x=list(x1,x2,x3,x4)
> model=lm(y~x1+x2+x3+x4)
> z=model$residuals
> check_periodicity(x,y,s)
Chi-squared test for detecting periodicity of coefficients






Chi-squared statistic: 68.60528  on df= 6 , p-value:  7.89701637415874e-13
> 
> ####""""""""""""""""""""""""" LSE method
> set.seed(4)
> s=2
> n=200
> m=n/s
> p=3
> mu=c(2,6)
> beta1=c(3,5)
> beta2=c(1,2.5)
> beta3=c(-2,1)
> x1=runif(n,0,5)
> x2=runif(n,0,10)
> x3=runif(n,0,15)
> y=rep(0,n)
> for (i in 1:s) {
+   q=seq(i,n,s)
+   y[q]=mu[i] + beta1[i] * x1[q] + beta2[i] * x2[q] + beta3[i] * x3[q] + rnorm(m,0,1)
+ }
> x=list(x1,x2,x3)
> lm_per(x,y,s)
Residuals:

-2.54963341465655-0.6268572651046620.02433322944958680.6396374150915482.920669384697835


Coefficients:

c(2.01783924645567, 5.75309920868897)c(3.07104872767226, 5.15277497228388)c(1.00521213892604, 2.47952546488017)c(-2.04168908550687, 0.98718191691937)c(0.963128009705219, 0.975882271186803)5


Root mean square error: 0.949937708348192
> 
> 
> ##""""""""""""""""""""""" AE
> set.seed(4)
> s=2
> n=200
> m=n/s
> p=3
> mu=c(2,6)
> beta1=c(3,5)
> beta2=c(1,2.5)
> beta3=c(-2,1)
> x1=runif(n,0,5)
> x2=runif(n,0,10)
> x3=runif(n,0,15)
> y=rep(0,n)
> for (i in 1:s) {
+   q=seq(i,n,s)
+   y[q]=mu[i] + beta1[i] * x1[q] + beta2[i] * x2[q] + beta3[i] * x3[q] + rsn(m,alpha = 10)
+ }
> x=list(x1,x2,x3)
> lm_per_AE(x,y,s)
Residuals:

-0.470631725170640.1400025684247150.4871983034858350.9210521813813412.650788651853465


Coefficients:

c(2.12500203467676, 6.24267402298002)c(2.96934921549602, 4.95748836239455)c(0.985143475664827, 2.51489863520411)c(-1.9698282085292, 0.99865589081242)c(0.467340237452444, 0.386634542672306)5


Root mean square error: 0.822466268155688
> 
> ######"""""""""""""""""""""""""" real data examples
> #####"
> #Data <- system.file("data", "weather_data_casablanca-settat.xlsx", package = "PerRegMod")
> #load("~/PerRegMod_4.4.1/PerRegMod/data/weather_data_casablanca_settat.RData")
> #Data=weather_data_casablanca_settat
> 
> #file_path <- system.file("data", "weather_data_casablanca_settat.xlsx", package = "PerRegMod")
> #Data <- readxl::read_excel(file_path)
> #file_path <- system.file("inst/extdata", "weather_data_casablanca_settat.xlsx", package = "PerRegMod")
> #Data <- readxl::read_excel(file_path)
> #Data=read_xlsx("inst/extdata/weather_data_casablanca_settat.xlsx")
> file_path <- system.file("extdata", "weather_data_casablanca_settat.xlsx", package = "PerRegMod")
> 
> # Check if the file path is found, then read the file
> if (file_path != "") {
+   Data <- readxl::read_excel(file_path)
+ } else {
+   stop("File not found. Please check that the file is correctly placed in 'inst/extdata' within the package directory.")
+ }
> months <- c('January','February','March','April','May','June',
+             'July','August','September','October','November','December')
> years <- c('2010','2011','2012','2013','2014','2015','2016',
+            '2017','2018','2019','2020')
> 
> #########""""""""""""for temperature
> y=rep(0,11 * 12)
> ly=length(years)
> lm=length(months)
> for (i in 1:ly ) {
+   for (j in 1:lm) {
+     y[(i - 1) * lm + j] = mean(Data$Temperature..C[Data$Year == years[i] &
+                                                      Data$Month == months[j] ]  )
+   }
+ 
+ }
> 
> #########""""""""""""for Humidity
> x1=rep(0,11 * 12)
> ly=length(years)
> lm=length(months)
> for (i in 1:ly ) {
+   for (j in 1:lm) {
+     x1[(i - 1) * lm + j]=mean(Data$Humidity..[Data$Year==years[i] &
+                                                 Data$Month==months[j] ]  )
+   }
+ 
+ }
> 
> #########""""""""""""for dew point
> x2=rep(0,11 * 12)
> ly=length(years)
> lm=length(months)
> for (i in 1:ly ) {
+   for (j in 1:lm) {
+     x2[(i - 1) * lm + j]=mean(Data$Dew.Point..C[Data$Year==years[i] &
+                                                   Data$Month==months[j] ]  )
+   }
+ 
+ }
> 
> #########""""""""""""for wind
> x3=rep(0,11 * 12)
> ly=length(years)
> lm=length(months)
> for (i in 1:ly ) {
+   for (j in 1:lm) {
+     x3[(i - 1) * lm + j]=mean(Data$Wind.Speed.km.h[Data$Year==years[i] &
+                                                      Data$Month==months[j] ]  )
+   }
+ 
+ }
> 
> ###""""""""""""""""""""""""example 1
> 
> #""""""""""check periodicity
> x=list(x2,x3)# dew point and wind
> n=length(y)
> s=12
> m=n/s
> check_periodicity(x,y,s)
Chi-squared test for detecting periodicity of coefficients






Chi-squared statistic: 139.35296  on df= 44 , p-value:  7.73558994637824e-12
> 
> #""""""""estimating parameters using LSE
> lm_per(x,y,s)
Residuals:

-3.79733393927142-0.725823019921449-0.02306340480145860.7315568062964723.387730659746615


Coefficients:

c(14.7685708035282, 14.4031208956541, 20.2313005468844, 16.6926932603836, 13.5494958194909, 5.18979739783482, 42.1295660598847, 37.8578774167528, 17.0023700763916, 32.701164138304, 26.9001054601609, 19.40800039521)c(0.365815587145837, 0.783058345168259, -0.0953134815151865, 0.893364714936179, 0.76741670948832, 1.60042623427687, -0.431044970184399, -0.247980083427971, 0.437179944161814, -0.306806324837901, -0.489266439126049, -0.0715023739455055)c(-0.344747045610824, -0.604913008116176, -0.210477052528003, -0.639931006866895, 0.00561041333596668, -0.304088120564935, -0.681611465595044, -0.453857332619767, 0.209452805975644, -0.503724563444341, -0.276029604005522, -0.278180551291081)c(0.77248044255742, 1.31027326304443, 1.01714214967979, 1.44246577642474, 2.07728957187128, 1.063017361327, 1.68553235454076, 1.31940302675612, 0.770235826903024, 1.72124395732565, 1.26231887692736, 1.15929220806617)5


Root mean square error: 1.15336557920339
> #""""""""estimating parameters using AE
> lm_per_AE(x,y,s)
Residuals:

5.044406433761259.189046572117612.744197585997415.103156991706728.96419656158395


Coefficients:

c(8.94063353330799, 3.63363558992864, 12.372573597533, 1.07415991565025, -1.74468550936581, -9.76069223473399, 17.2840835056552, 24.3599865827471, 9.0579115438521, 19.1111085363004, 15.0133284842706, 6.77877760775152)c(0.367149246381972, 0.818821687854316, -0.0439407566145752, 0.808259636423481, 0.658553506058031, 1.65750273571941, -0.495887852154496, -0.200481265982059, 0.447995794574028, -0.27559827428159, -0.486378638205221, -0.0716520396070216)c(-0.34617813408455, -0.627715985470845, -0.233394686282741, -0.586823779066773, 0.0816155287819667, -0.3424465281827, -0.672041515247122, -0.473600572375121, 0.202956994486457, -0.570506060304705, -0.2823016970107, -0.276095759896722)c(-8.4356823262697, -46.3727922652315, -20.1689379782247, -82.1568714075151, -166.695603745412, -42.7829428275578, -180.176702765697, -59.3801123165356, -11.6529914143607, -101.279643867707, -47.5739305852113, -42.7480310305977)5


Root mean square error: 13.9178614908708
> 
> ####"""""""""""""""""""""Example 2
> ##""""""check periodicity
> x=list(x1,x3)# humidity and wind
> n=length(y)
> s=12
> m=n/s
> check_periodicity(x,y,s)
Chi-squared test for detecting periodicity of coefficients






Chi-squared statistic: 1131.78449  on df= 44 , p-value:  0
> 
> #""""""""estimating parameters using LSE
> lm_per(x,y,s)
Residuals:

-2.36837425068774-0.47925980656605-0.08392215182808370.5443225929583383.561057739334325


Coefficients:

c(14.9795538851734, 20.709623357528, 29.0594613482795, 46.4671758299881, 43.1574955641698, 50.8715585678433, 48.4717948093222, 45.1517388982596, 44.8284606355821, 42.2374952967589, 32.4659909188156, 24.9660007299092)c(0.0333039895713665, -0.089700727361186, -0.156964861935598, -0.35488369176146, -0.306016266750712, -0.378169315200616, -0.322813820225341, -0.250711745020876, -0.262107619899439, -0.232873014851319, -0.163448047967214, -0.106497815158633)c(-0.313530753205975, 0.122883069056731, -0.0611943472412457, -0.268506644104732, -0.115093015878652, -0.225891886724934, -0.140465229259265, -0.179925756961408, -0.236498466855422, -0.458102115385911, -0.290638068680655, -0.155003997874916)c(0.930410693291833, 1.9082697146199, 0.644805992172765, 0.825675541766766, 1.51030660365062, 0.999130413104471, 1.00969797390427, 0.6856030837865, 0.573713802598791, 0.635772527297887, 0.752709916254286, 0.962006933520253)5


Root mean square error: 0.873824239594352
> #""""""""estimating parameters using AE
> lm_per_AE(x,y,s)
Residuals:

-75.2219473511798-26.4136489314989-11.98180628193921.2174302147825617.24558034226645


Coefficients:

c(14.8193117078783, 21.0334794465951, 29.0848244659445, 46.3464052686244, 42.908497411699, 51.0071886196358, 48.5944392394843, 45.1792024437559, 44.8726464216715, 42.2474939506746, 32.3949469786859, 24.9507886749815)c(0.361552537789099, 0.967482485132604, -0.418684133490646, -0.352239617500968, -0.528517014710358, 0.620543134141942, 0.108567724235891, 0.0166750255750096, -0.31564061087719, -0.140318641182971, 0.0368061072744772, 0.539494451284724)c(-0.656014155087793, -0.962244144231961, 0.183079546831297, -0.276514469957811, 0.152993658686196, -1.18994235840608, -0.552599028233534, -0.433678035139622, -0.180575018155406, -0.550130403171244, -0.497055454095078, -0.827190943470303)c(1.01933768169097, 2.05471318871494, 0.395864311292907, 0.778368688596071, 3.10222591411073, 0.872479638463224, 0.891562058774985, 0.462078939625592, 0.312672781007426, 0.399335843897158, 0.611414912679856, 0.911178703252318)5


Root mean square error: 29.2058297009007
> 
> 
> proc.time()
   user  system elapsed 
   6.20    0.40    6.59