> library(tseries)
> library(Rssa)

> hujan=read.csv(file.choose())
> Curah_hujan<-ts(hujan[ , ,12], start=c(2014,1),frequency=12)
> acf(Curah_hujan)
> pacf(Curah_hujan)
> plot(Curah_hujan)
> adf.test(Curah_hujan)

> inhujan=read.csv(file.choose())
> data=as.matrix(inhujan)
> n<-nrow(data)
> L=36
> K=n-L+1
> lintasan<-matrix(,nrow=L,ncol=K)
> for(i in 1:K){
+ if(i==1){
+ lintasan[,i]<-data[i:L]
+ }else{
+ lintasan[,i]<-data[i:(L+(i-1))]}}

> S_simetris=lintasan%*%(t(lintasan))
> e=eigen(S_simetris)
> eigenvalue=as.matrix(e$values)
> plot(eigenvalue,type="l",col="blue")

> singular_value=sqrt(eigenvalue)
> plot(singular_value,type="l",col="blue")

> eigenvector=as.matrix(e$vectors)
> plot(eigenvector[,1],eigenvector[,2],type="l",col="blue")

> evecs=1/singular_value
> Vbaru<-list()
> for(i in 1:L){
+ Vbaru[[i]]<-evecs[i]*((t(lintasan))%*%((eigenvector[,i])))}
> Vi<-do.call(cbind,Vbaru)

> opar=par(no.readonly=TRUE)
> par(mfrow=c(3,3))
> plot(Vi[,1],type="l",col="blue")

> s<-ssa(hujan,L=36)
> plot(wcor(s))

> Xbaru<-list()
> for(i in 1:L){
+ Xbaru[[i]]<-singular_value[i]*((eigenvector[,i]))%*%as.matrix(t(Vi[,i]))}
> X1=(Xbaru[[1]])
> X2=Reduce('+',Xbaru[c(2,4)])
> X3=(Xbaru[[3]])
> X4=(Xbaru[[5]])
> X5=(Xbaru[[6]])
> X6=(Xbaru[[7]])
> X7=(Xbaru[[8]])
> X8=(Xbaru[[9]])
> X9=(Xbaru[[10]])
> X10=(Xbaru[[11]])
> X11=(Xbaru[[12]])
> X12=(Xbaru[[13]])
> X13=(Xbaru[[14]])
> X14=(Xbaru[[15]])
> X15=(Xbaru[[16]])
> X16=(Xbaru[[17]])
> X17=(Xbaru[[18]])
> X18=(Xbaru[[19]])
> X19=(Xbaru[[20]])
> X20=(Xbaru[[21]])
> X21=Reduce('+',Xbaru[c(22,23)])
> X22=Reduce('+',Xbaru[c(24,25)])
> X23=(Xbaru[[26]])
> X24=(Xbaru[[27]])
> X25=(Xbaru[[28]])
> X26=(Xbaru[[29]])
> X27=(Xbaru[[30]])
> X28=(Xbaru[[31]])
> X29=(Xbaru[[32]])
> X30=(Xbaru[[33]])
> X31=Reduce('+',Xbaru[c(34,35)])
> X32=(Xbaru[[36]])

> recon1<-reconstruct(s,groups=list(1,c(2,4),3,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,c(22,23),c(24,25),26,27,28,29,30,31,32,33,c(34,35),36))
> f<-rforecast(s,groups=list((1:32)),only.new=FALSE,len=36)

> outhujan=read.csv(file.choose())
> aktual=as.matrix(outhujan)
> ramal=f[85:96]
> par(mfrow=c(1,1))
> plot(aktual, type="o", col="blue",lwd=2,xlab="month",ylab="Curah Hujan (mm)", main="Data Aktual Dan Ramalan")
> lines(ramal,col="red",lwd=2)

> PE=array(NA,dim=c(12))
> for(i in 1:12){
+ PE[i]=abs((aktual[i]-ramal[i])/aktual[i])}
> PE
> MAPE=mean(PE, na.rm=TRUE)
> MAPE