R Under development (unstable) (2024-11-11 r87319 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")
Installing package into 'D:/temp/2024_11_11_23_35_02_18960/RtmpeEwbTC/RLIBS_272a840be63f6'
(as 'lib' is unspecified)
trying URL 'https://cran.r-project.org/bin/windows/contrib/4.5/expm_1.0-0.zip'
Content type 'application/zip' length 222460 bytes (217 KB)
==================================================
downloaded 217 KB

package 'expm' successfully unpacked and MD5 sums checked

The downloaded binary packages are in
	D:\temp\2024_11_11_23_35_02_18960\RtmpG4Ataw\downloaded_packages
> install.packages("sn", repos = "https://cran.r-project.org")
Installing package into 'D:/temp/2024_11_11_23_35_02_18960/RtmpeEwbTC/RLIBS_272a840be63f6'
(as 'lib' is unspecified)
trying URL 'https://cran.r-project.org/bin/windows/contrib/4.5/sn_2.1.1.zip'
Content type 'application/zip' length 1845934 bytes (1.8 MB)
==================================================
downloaded 1.8 MB

package 'sn' successfully unpacked and MD5 sums checked

The downloaded binary packages are in
	D:\temp\2024_11_11_23_35_02_18960\RtmpG4Ataw\downloaded_packages
> install.packages("readxl", repos = "https://cran.r-project.org")
Installing package into 'D:/temp/2024_11_11_23_35_02_18960/RtmpeEwbTC/RLIBS_272a840be63f6'
(as 'lib' is unspecified)
trying URL 'https://cran.r-project.org/bin/windows/contrib/4.5/readxl_1.4.3.zip'
Content type 'application/zip' length 1203600 bytes (1.1 MB)
==================================================
downloaded 1.1 MB

package 'readxl' successfully unpacked and MD5 sums checked

The downloaded binary packages are in
	D:\temp\2024_11_11_23_35_02_18960\RtmpG4Ataw\downloaded_packages
> 
> #install.packages("expm")
> #install.packages("sn")
> #install.packages("readxl")
> 
> 
> 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: 48.56717  on df= 6 , p-value:  9.10120778740975e-09> 
> ####""""""""""""""""""""""""" 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:
       0%       25%       50%       75%      100% 
-2.549633 -0.626857  0.024333  0.639637  2.920669 

Coefficients:
     intercept x...1... X...2... X...3...      sd
s= 1    2.0178   3.0710   1.0052 -2.04169 0.96313
s= 2    5.7531   5.1528   2.4795  0.98718 0.97588

Root mean square error:  0.9499377> 
> 
> ##""""""""""""""""""""""" 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%      25%      50%      75%     100% 
-0.47063  0.14000  0.48720  0.92105  2.65079 

Coefficients:
     intercept x...1... X...2... X...3...      SD
s= 1    2.1250   2.9693  0.98514 -1.96983 0.46734
s= 2    6.2427   4.9575  2.51490  0.99866 0.38663

Root mean square error:  0.8224663> 
> ######"""""""""""""""""""""""""" 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:
       0%       25%       50%       75%      100% 
-3.797334 -0.725823 -0.023063  0.731557  3.387731 

Coefficients:
      intercept  x...1...   X...2...      sd
s= 1    14.7686  0.365816 -0.3447470 0.77248
s= 2    14.4031  0.783058 -0.6049130 1.31027
s= 3    20.2313 -0.095313 -0.2104771 1.01714
s= 4    16.6927  0.893365 -0.6399310 1.44247
s= 5    13.5495  0.767417  0.0056104 2.07729
s= 6     5.1898  1.600426 -0.3040881 1.06302
s= 7    42.1296 -0.431045 -0.6816115 1.68553
s= 8    37.8579 -0.247980 -0.4538573 1.31940
s= 9    17.0024  0.437180  0.2094528 0.77024
s= 10   32.7012 -0.306806 -0.5037246 1.72124
s= 11   26.9001 -0.489266 -0.2760296 1.26232
s= 12   19.4080 -0.071502 -0.2781806 1.15929

Root mean square error:  1.153366> #""""""""estimating parameters using AE
> lm_per_AE(x,y,s)
Residuals:
     0%     25%     50%     75%    100% 
 5.0444  9.1890 12.7442 15.1032 28.9642 

Coefficients:
      intercept  x...1...  X...2...        SD
s= 1     8.9406  0.367149 -0.346178   -8.4357
s= 2     3.6336  0.818822 -0.627716  -46.3728
s= 3    12.3726 -0.043941 -0.233395  -20.1689
s= 4     1.0742  0.808260 -0.586824  -82.1569
s= 5    -1.7447  0.658554  0.081616 -166.6956
s= 6    -9.7607  1.657503 -0.342447  -42.7829
s= 7    17.2841 -0.495888 -0.672042 -180.1767
s= 8    24.3600 -0.200481 -0.473601  -59.3801
s= 9     9.0579  0.447996  0.202957  -11.6530
s= 10   19.1111 -0.275598 -0.570506 -101.2796
s= 11   15.0133 -0.486379 -0.282302  -47.5739
s= 12    6.7788 -0.071652 -0.276096  -42.7480

Root mean square error:  13.91786> 
> ####"""""""""""""""""""""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:
       0%       25%       50%       75%      100% 
-2.368374 -0.479260 -0.083922  0.544323  3.561058 

Coefficients:
      intercept  x...1...  X...2...      sd
s= 1     14.980  0.033304 -0.313531 0.93041
s= 2     20.710 -0.089701  0.122883 1.90827
s= 3     29.059 -0.156965 -0.061194 0.64481
s= 4     46.467 -0.354884 -0.268507 0.82568
s= 5     43.157 -0.306016 -0.115093 1.51031
s= 6     50.872 -0.378169 -0.225892 0.99913
s= 7     48.472 -0.322814 -0.140465 1.00970
s= 8     45.152 -0.250712 -0.179926 0.68560
s= 9     44.828 -0.262108 -0.236498 0.57371
s= 10    42.237 -0.232873 -0.458102 0.63577
s= 11    32.466 -0.163448 -0.290638 0.75271
s= 12    24.966 -0.106498 -0.155004 0.96201

Root mean square error:  0.8738242> #""""""""estimating parameters using AE
> lm_per_AE(x,y,s)
Residuals:
      0%      25%      50%      75%     100% 
-75.2219 -26.4136 -11.9818   1.2174  17.2456 

Coefficients:
      intercept  x...1... X...2...      SD
s= 1     14.819  0.361553 -0.65601 1.01934
s= 2     21.033  0.967482 -0.96224 2.05471
s= 3     29.085 -0.418684  0.18308 0.39586
s= 4     46.346 -0.352240 -0.27651 0.77837
s= 5     42.908 -0.528517  0.15299 3.10223
s= 6     51.007  0.620543 -1.18994 0.87248
s= 7     48.594  0.108568 -0.55260 0.89156
s= 8     45.179  0.016675 -0.43368 0.46208
s= 9     44.873 -0.315641 -0.18058 0.31267
s= 10    42.247 -0.140319 -0.55013 0.39934
s= 11    32.395  0.036806 -0.49706 0.61141
s= 12    24.951  0.539494 -0.82719 0.91118

Root mean square error:  29.20583> 
> 
> proc.time()
   user  system elapsed 
   7.32    0.48    8.73