##########################################
##
## calculate forecast of covid-19 
## using prophet
##
## 14.4.2022
## v 0000.0006
##
##

#
#install.packages("rvest")
#install.packages("readtext")
#install.packages("stringi")
#install.packages("datamart")
#install.packages("XML")
#install.packages("svglite")
#install.packages("ggplot2")
#install.packages("tidyr")
#install.packages("stringr")
#install.packages("stringi")
#install.packages("tibble")

Sys.setlocale("LC_ALL","Finnish")
options(encoding = "UTF-8")

library(ggplot2)
library(svglite)

library(rvest)
library(readtext)
library(stringi)
library(stringr)
library(datamart)
library(XML)
library(jsonlite)
library(rjstat)

library(tibble)
library(caTools)
library(mgcv)
library(repmis)
library(lubridate)
library(tidyverse)
library(tidyr)
library(dplyr)
library(data.table)

#library(covid19.analytics)
library(R0)
library(EpiEstim)
library(prophet)

## choices
## 1 finnish wiki data, 2 aggregated cases data 
## 3 solanpaa finnish data 4 thl cube json data
#

load_data_from=3


yala=0
yyla=15000


beginday1='01/11/2021'

forecastendday1<-"2022/07/01"

## 1 from finnish wiki, 2 cases from net, 3 from net 2
plottaa=1  ## must be 1
tulosta_svg=1     # plot to out svg 0, 1 of 2

tulosfilee1="/Users/himot/akor1/prophet1.svg"

## ggolot smooth curves pars
spanni=0.2
#spanni=0.5

metodi="loess"

## NOTE date limits change this

#datelimits1=c('1/3/2020', '9/11/2020')
## display date limits

today=Sys.Date()-12

#print(today)

today1=format(today, "%d/%m/%Y")  
today2=format(today, "%Y/%m/%d")  

print(today1)
#print(today2)

#stop(-1)

datelimits1=c(beginday1, today1)

print(today1)
print(datelimits1)

#stop(-1)


## dates of dataset

paivat1=seq(as.Date("2020/4/1"), as.Date(today2), "days")



calculate_r0 <- function(time1, time2, val1, val2) 
{
	td=time2-time1
	gr0<-log(val2/val1)
	gr=gr0/td
	td = log(2)/gr
	tau<-5.0
	k<-log(2.0)/td
	r0<-exp(k*tau)
	return(r0)
}


moving_average <- function(x, w, FUN, ...) 
{

  if (w < 1) {
    stop("Window length: mustbe greater than 0")
  }
  output <- x
  for (i in 1:length(x)) {
    lower_bound <- i - w + 1
    if (lower_bound < 1) {
      output[i] <- NA_real_
      ## !!! assume NA 0
      output[i] <- 0
    } else {
      output[i] <- FUN(x[lower_bound:i, ...])
    }
  }
  return (output)
}

calculate_multiple_r0 <- function(daata1) {

  lenu1<-length(daata1)
   
  daata2<-1:lenu1
  
	for (n in 2:lenu1){
		valju1=daata1[n-1]
		valju2=daata1[n]
		timex1=0
		timex2=1

		r0<-calculate_r0(0, 1, valju1, valju2)
		daata2[n]<-r0

		#print (r0)
		}
	return(daata2)
	
}


load_data_from_finnish_wiki<-function()
{

url1="https://fi.wikipedia.org/wiki/Suomen_koronaviruspandemian_aikajana"
destfile1="./ward0.txt"

download.file(url1, destfile1)
texti000<-readtext(destfile1) 
texti0<-texti000$text

etsittava1="1. huhtikuuta 2020 alkaen"
len1=nchar(texti0)
k1=regexpr(pattern=etsittava1, texti0)
k1b=len1-k1
texti1=strtail(texti0,k1b)
sink("out1.txt")
print (texti1)
sink()
etsittava2="</table>"
k2=regexpr(pattern=etsittava2, texti1)
texti2=strhead(texti1,k2)

sample1<-minimal_html(texti2)
tabu1 <- html_table(sample1, fill=TRUE)[[1]]
colnames(tabu1) <- c("V1","V2", "V3","V4", "V5","V6", "V7","V8" )

#print(tabu1)

sairaalassa00<-tabu1$V4
sairaalassa=as.integer(sairaalassa00)

teholla00<-tabu1$V5
teholla=as.integer(teholla00)

uusiatapauksia00<-tabu1$V3
uusiatapauksia0<-gsub(" ", "", uusiatapauksia00)
uusia_tapauksia=as.integer(uusiatapauksia0)

uusiakuolleita00<-tabu1$V7
uusiakuolleita1=as.integer(uusiakuolleita00)

uusiakuolleita2<-uusiakuolleita1
uusiakuolleita2[uusiakuolleita2<0]<-0
uusia_kuolleita<-uusiakuolleita2

toipuneita00<-tabu1$V8
toipuneita01<-gsub(" ", "", toipuneita00)
toipuneita0<-gsub("[^0-9.-]", "", toipuneita01)
toipuneita=as.integer(toipuneita0)

tapauksia00<-tabu1$V2
tapauksia01<-gsub(" ", "", tapauksia00)
tapauksia0<-gsub("[^0-9.-]", "", tapauksia01)

tapauksia=as.integer(tapauksia0)

kuolleita00<-tabu1$V6
kuolleita=as.integer(kuolleita00)


aktiivisia_tapauksia=tapauksia-kuolleita-toipuneita


#print (paivat1)
#print (teholla)
#print (sairaalassa)
#print (tapauksia)
#print (kuolleita)
#print (toipuneita)
#print (uusia_tapauksia)
#print (uusia_kuolleita)
#plot(paivat1,aktiivisia_tapauksia)

#xy<-data.frame(paivat1, sairaalassa)

xy<-data.frame(paivat1, uusia_tapauksia)


xyz<-data.frame(paivat1, sairaalassa, teholla) 
dfout1<-data.frame(paivat1, aktiivisia_tapauksia, uusia_tapauksia, sairaalassa, teholla, uusia_kuolleita ) 

names(dfout1)<-c("Pvm", "Aktiivisia_tapauksia","Uusia_tapauksia", "Sairaalassa", "Teholla", "Uusia_kuolleita")

write.csv2(dfout1, "./sairaalassa.csv",row.names=FALSE )

return(xy)
}


load_data_from_aggregated<-function()
{
###################################################################################
###################################################################################
## fetch the data

dfine <- read.csv(file = 'https://datahub.io/core/covid-19/r/countries-aggregated.csv')

#head(dfine)
#class(dfine)

#tail(dfine, 5)

dfinland <- dfine[ which(dfine$Country=='Finland'), ]

#head(dfinland)

kols <- c("Date", "Confirmed","Recovered","Deaths")

tapaukset <- dfinland[kols]

#head(tapaukset)

len1=nrow(tapaukset)

#len1

len2=len1-1

len3=len2

confirmed<-tapaukset$Confirmed
deaths<-tapaukset$Deaths

dailycases <- vector()
dailycases <- c(dailycases, 0:(len2))
dailydeaths <- vector()
dailydeaths <- c(dailydeaths, 0:(len2))

m=0
dailycases[1]<-tapaukset$Confirmed[1]
#dailydeaths[1]<-tapaukset$Deaths[1]
dailydeaths[1]<-0

#confirmed
#deaths

m=1
for(n in 2:(len3+1)) {
	
	a<-confirmed[n]
	b<-confirmed[m]
	#print (a)
	#print (b)
	cee<- (a-b)
	#print(cee)
	dailycases[n]=cee	
	m=m+1
}

mm=1
for(nn in 2:(len3+1)) {
	
	aa<-deaths[nn]
	bb<-deaths[mm]
	#print ("_")
	#print (aa)
	#print (bb)
	ceb=aa-bb
	#if (ceb<0) ceb=0
	#print(ceb)
	dailydeaths[nn]=ceb	
	mm=mm+1
}

#deaths

#dailycases

#dailydeaths

dfout1<-dfinland

#print(nrow(dfinland))
#print(length(dailydeaths))

dfout1 <- cbind(dfout1, data.frame(dailycases))
dfout1 <- cbind(dfout1, data.frame(dailydeaths))

#head(dfout1)

dfout2<-within(dfout1, rm(Country))

names(dfout2) <- c('Date','Confirmed','Recovered','Deaths', 'DailyConfirmed','DailyDeaths')

#head(dfout2)

write.csv2(dfout2, "/Users/himot/akor1/finland_data1.csv");

daate1<-dfout2$Date
dailydeaths1<-dfout2$DailyDeaths
dailycases1<-dailycases

#daate1

#daate2<-gsub("2020-", "", daate1)
daate2<-daate1

leenu<-length(daate2)

#alkupvm<-50
alkupvm<-1

daate3<-daate2[alkupvm:leenu]
dailydeaths3<-dailydeaths1[alkupvm:leenu]
dailycases3<-dailycases1[alkupvm:leenu]

#daate3
#dailydeaths3

    
 # barplot(dailydeaths3, main="Koronaviruskuolemat päivittäin vuonna 2020",
 # names.arg=daate3) 
  
      
  dataf1 <- data.frame("Date" = daate3, "Paivitt_kuolemat"=dailydeaths3)
# str(dataf1)    
 
   dataf2 <- data.frame("Date" = daate3, "Paivitt_tapaukset"=dailycases3)
# str(dataf2)  
 
 write.csv(dataf1, "/Users/himot/akor1/dailydeaths1.csv", row.names=T)
 write.csv(dataf2, "/Users/himot/akor1/dailycases1.csv", row.names=T)

 indf1 <- read.csv(file = '/Users/himot/akor1/dailycases1.csv')
 #head(indf1)

 cases1<-indf1$Paivitt_tapaukset
 dates1<-indf1$Date

 len1=length(cases1)

 dates2<-as.Date(dates1)
 paivat<-1:len1	
 
 xy<-data.frame(daate3, dailycases3)
 
	
}




calculate_r0_with_r0<-function(xy2)
{
	
	## calculate r0 w/r0 package
	dates<-as.Date(xy2$Dates)
	cases<-as.integer(xy2$Cases)
	
	cases[is.na(cases)] <- 1
	cases[(cases<0)] <- cases*-1
	cases[cases==0] <- 1
	
	nummeros<-1:length(dates)
	num<-cases
	#names<-nummeros
	names<-dates
	lenu=length(dates)
	
	bekini=as.Date(dates[1])
	enti=as.Date(dates[lenu])
	
	#print(bekini)
	#print(enti)
	
	#stop(-1)
	
	#enti=lenu
		
	#bekini=enti*0+1
	
	#enti=as.integer(enti)
	#bekini=as.integer(bekini)
	
	
	df1 <- setNames(num, names)
	
	mGT<-generation.time("gamma", c(3, 1.5))
	#TD <- est.R0.TD(df1, mGT, begin=1, end=length(dates), nsim=200)
	#TD <- est.R0.TD(df1, mGT, begin=bekini, end=enti, nsim=200)
	TD <- est.R0.TD(df1, mGT, begin=bekini, end=enti, nsim=200)
	TD.5D <- smooth.Rt(TD, 5)
	paivat1<-TD.5D$epid$t 
	paivat2<-as.Date(paivat1)
	r0t1<-TD.5D$R
	conf1<-TD.5D$conf.int
	
	xypaluu<-data.frame(paivat1,r0t1)
	names(xypaluu)<-c("paivat","r0")
	return(xypaluu)
}


calculate_r0_with_epiestim<-function(xy2)
{
	
	## calculate r0 w/r0 package
	dates<-as.Date(xy2$Dates)
	cases<-as.integer(xy2$Cases)
	nummeros<-1:length(dates)
	num<-cases
	#names<-nummeros
	names<-dates
	lenu=length(dates)
	
	cases[is.na(cases)] <- 1
	cases[(cases<0)] <- cases*-1
	cases[cases==0] <- 1
	
	incid<-cases
	
	bekini=as.Date(dates[1])
	enti=as.Date(dates[lenu])
	
	config<-make_config( list(mean_si = 2.6,std_si = 1.5) )

	res<-estimate_R(incid,method="parametric_si", config = config)   
                                 
	#plot(res)
	resr<-res$R

	str(resr)

	meanr<-resr$Mean 
	medianr<-resr$Median 
	quantile95<-resr$Quantile.0.95
	quantile05<-resr$Quantile.0.05
	quantile75<-resr$Quantile.0.75
	quantile25<-resr$Quantile.0.25
	meanr

	daydexes<-resr$t_start 

	daydexes
	
	
	#plot(daydexes, meanr)

	dayss<-as.Date(dates[daydexes])

	print (dayss)
	#stop(-1)
	#plot(dayss, meanr)
	
	xypaluu<-data.frame(dayss,meanr)
	names(xypaluu)<-c("paivat","r0")
	return(xypaluu)
}


calculate_r0_with_simple_exponent_moving_average<-function(xy2, madays1, madays2)
{
	
	## calculate r0 w/r0 package
	dates<-as.Date(xy2$Dates)
	cases<-as.integer(xy2$Cases)
	nummeros<-1:length(dates)
	num<-cases
	#names<-nummeros
	names<-dates
	lenu=length(dates)
	
	cases[is.na(cases)] <- 1
	cases[(cases<0)] <- cases*-1
	cases[cases==0] <- 1
 
	# compute a MA(7)
	ma1<-moving_average(cases,madays1,mean)
	r0t1<-calculate_multiple_r0(ma1)
	r0avg1<-moving_average(r0t1, madays2, mean)
	xypaluu<-data.frame(dates,r0t1)
	
	#plot(r0t1)
	#print (r0t1)
	#stop(-1)
	
	names(xypaluu)<-c("paivat","r0")
	
	return(xypaluu)
	
}


forecast_profet<-function(xy2, futuredays)
{
		## calculate r0 w/r0 package
	ds<-as.Date(xy2$Dates)
	y<-as.integer(xy2$Cases)
	nummeros<-1:length(ds)

	lenu=length(ds)
	
	
	
	df<-data.frame(ds,y)


	m <- prophet(df)
	
	future <- make_future_dataframe(m, periods = futuredays)
	
	forecast <- predict(m, future)


    #str(future)
    #str(forecast)
    
    futu_days=future$ds
    futu_trendi=forecast$trend
    futu_trendi_upper=forecast$trend_upper
    futu_trendi_lower=forecast$trend_lower
    futu_yhat=forecast$yhat
    futu_yhat_upper=forecast$yhat_upper
    futu_yhat_lower=forecast$yhat_lower
    futu_weekly=forecast$weekly
    futu_weekly_upper=forecast$weekly_upper
    futu_weekly_lower=forecast$weekly_lower
    
    xypaluu<-data.frame(as.Date(futu_days),futu_yhat)
	# xypaluu<-data.frame(as.Date(futu_days),futu_weekly)
	#plot(r0t1)
	#print (r0t1)
	#stop(-1)
	
	names(xypaluu)<-c("paivat","r0")
	return(xypaluu)
    
}


lataa_thl_tapaukset_kuolleet<-function()
{
    url1<-"https://sampo.thl.fi/pivot/prod/fi/epirapo/covid19case/fact_epirapo_covid19case.json?row=measure-492118&column=dateweek20200101-508804L"
	cube1 <- fromJSONstat(url1, naming = "label", use_factors = F, silent = T)
	res01 <- cube1[[1]]
	#res00
	url2<-"https://sampo.thl.fi/pivot/prod/fi/epirapo/covid19case/fact_epirapo_covid19case.json?row=measure-444833&column=dateweek20200101-508804L"
	cube2 <- fromJSONstat(url2, naming = "label", use_factors = F, silent = T)
	res02 <- cube2[[1]]
	#res02

	#stop (-1)
	paiva=as.Date(res01$dateweek20200101)
	kuolleet=as.integer(res01$value)
	tapaukset=as.integer(res02$value)

	kuolin_prosentit=kuolleet/tapaukset
	kuolin_prosentit=kuolin_prosentit*10000
	kuolin_prosentit=as.integer(kuolin_prosentit)
	kuolin_prosentit=as.double(kuolin_prosentit)
	kuolin_prosentit=kuolin_prosentit/100.0

	#print (paiva)
	#print (kuolleet)
	#stop(-1)
	#print (tapaukset)
	#print (kuolin_prosentit )

	df1<-data.frame(paiva,tapaukset, kuolleet, kuolin_prosentit) 

	names(df1)<-c("Paiva", "Tapauksia", "Kuolleita", "Kuolinprosentti")
	#write.csv2(df1, "./kuolleet_ikaryhmittain.csv", sep = ";" )

	write.csv(df1, "./thl_tapaukset_kuolleet.csv")
	
	xy0<-data.frame(paiva, tapaukset)
	names(xy0)<-c("Dates", "Cases")
	xy<-na.omit(xy0)
	
	#return(df1)
	
	
	}
  

download_solanpaa_finnish_data<-function()
{
	solanpaa_fi="https://covid19.solanpaa.fi/data/fin_cases.json"
	cache_file="solanpaa_fi.json"

	download.file(solanpaa_fi, cache_file)

	j1 <- fromJSON(cache_file)

  ## maybe errori
	dates<-as.Date(j1$date)
	
	dailycases<-j1$new_cases
	dailydeaths<-j1$new_deaths

    dataf1 <- data.frame("Date" = dates, "Paivitt_kuolemat"=dailydeaths) 
    dataf2 <- data.frame("Date" = dates, "Paivitt_tapaukset"=dailycases)
 
	write.csv(dataf1, "./dailydeaths1.csv", row.names=T)
	write.csv(dataf2, "./dailycases1.csv", row.names=T)
	
	xy0<-data.frame(dates, dailycases)
	names(xy0)<-c("Dates", "Cases")
	xy<-na.omit(xy0)
	
	return(xy)	
	
}	




if(load_data_from==1)
{
	xy<-load_data_from_finnish_wiki()
	print (xy)
}

if(load_data_from==2)
{
	xy<-load_data_from_aggregated()
}

if(load_data_from==3)
{
	xy<-download_solanpaa_finnish_data()
}

if(load_data_from==4)
{
	xy<-lataa_thl_tapaukset_kuolleet()
}





 names(xy)<-c("Dates","Cases")


 #print (xy)

 #print (beginday1)
 select_datelimit_begin=as.Date(beginday1,format="%d/%m/%Y")
 select_datelimit_end=as.Date(today1, format= "%d/%m/%Y")

 #format(select_datelimit_begin, "%Y-%m-%d")

 #print(select_datelimit_begin)

 #print(select_datelimit_end)
 #2020-12-16

 #xy2<-xy[xy$Dates >= select_datelimit_begin & xy$Dates <= select_datelimit_end]
  
  #xy20<-xy[xy$Dates >= select_datelimit_begin,]
  #xy2<-xy20[xy2$Dates <= select_datelimit_end, ] 
  #str(xy2)
#  xy2 <- subset(xy, Dates > "2015-01-01" & Dates < "2015-12-31")
##  xy2 <- subset(xy, Dates > select_datelimit_begin & Dates < select_datelimit_end )

	xyt<-xy
	setDT(xyt)  

	#testset[date %between% c("2013-08-02", "2013-11-01")]
	#xyt2<-xyt[Dates %between% c(select_datelimit_begin,select_datelimit_end)]
	#xyt2 <- subset(xyt, Dates > "2021-07-01" & Dates < "2021-12-04")
	xyt20 <- subset(xyt, Dates > select_datelimit_begin )
	xyt2  <- subset(xyt20, Dates < select_datelimit_end )
#	xy2<-xy[xy$Dates >= select_datelimit_begin,]
 
	print(select_datelimit_begin)
	print(select_datelimit_end)
	print(xyt2)
	
	xy2<-as.data.frame(xyt2)
 
#  stop(-1)
  
 cases1<-xy2$Cases
 dates1<-xy2$Dates
 
 xy3<-data.frame( as.Date(dates1),as.integer(cases1) )
 names(xy3)<-c("Dates", "Cases")

 len1=length(cases1)

 dates2<-as.Date(dates1)
 paivat<-1:len1

 
 ## test code
 arrat0<-calculate_r0_with_simple_exponent_moving_average(xy2, 14,7)
 
 arrat1<-calculate_r0_with_r0(xy2)
 arrat2<-calculate_r0_with_epiestim(xy2)
 
  #plot(arrat$paivat, arrat$r0)

  arrat<-arrat2
  

 
 
 #str(arrat)
 #head(arrat)


  sarrat1<-arrat1
  sarrat2<-sarrat1
  names(sarrat1)<-c("Dates","Cases")
   
  
  datelimits2=c(today1, as.Date(forecastendday1,"%Y/%m/%d"))
  datelimits3=c(as.Date(beginday1, "%d/%m/%Y" ), as.Date(forecastendday1,"%Y/%m/%d"))

  daysek1<-seq(today, as.Date(forecastendday1), "days")
  lendaysek1<-length(daysek1)

 bed1<-as.Date(beginday1, "%d/%m/%Y")
 
 #end1<-as.Date(forecastendday1)+1
   end1<-as.Date(forecastendday1)-1

  print(bed1)
  print (end1)

 # stop(-1)
  

  daysek2<-seq(bed1, end1, "days")
  daysek3<-seq(as.Date(beginday1,"%d/%m/%Y" ),as.Date(forecastendday1), "days")




  ###
  ###
  basedata<-xy2
  names(basedata)<-c("ds", "y")
  #farrat1<-forecast_profet(xy2, lendaysek1)
  
   profeta <- prophet(basedata,interval.width = 0.1, yearly.seasonality =0, weekly.seasonality = TRUE)
   future <- make_future_dataframe(profeta, periods = lendaysek1)
   forecast1 <- predict(profeta, future)
    
    plot(profeta, forecast1)
    
    
  #  print(forecast1)
    
    yhat1<-forecast1$yhat
        
    #lower1<-yhat1-forecast1$trend_lower/2
    #upper1<-yhat1+forecast1$trend_upper/2
    
   lower1<-yhat1-forecast1$yhat_lower
    upper1<-yhat1+forecast1$yhat_upper
    
  #  tail(forecast1[c(‘ds’, ‘yhat’, ‘yhat_lower’, ‘yhat_upper’)])
  # head(forecast1)
   
  #  print(daysek2)
    
  #  stop(-1)
    
    	print(daysek2)
  #	print(yhat1)
  #		print(lower1)
 # 	print(upper1)
  #	print(length(daysek2))
  #	print(length(yhat1))
  	
  
    farrat1<-data.frame(daysek2,yhat1,lower1, upper1)
    
    print(farrat1)
 
 # stop(-1)
    
    
  names(farrat1)<-c("Dates", "Cases", "Lower", "Upper")
  
  print(datelimits3)
  
  
  #print(farrat1)
  
   
  
  
 # stop(-1)
  
  
  #plot(arrat$paivat, arrat$r0)
 # plot(farrat1$paivat, farrat1$r0)
 # lines(farrat1$paivat, farrat1$r0  )
 # lines(sarrat2$paivat, sarrat1$r0  )
 
#	stop(-1)
 
 arrat<-farrat1
 
 
 
  names(arrat)<-c("Dates","Forecast")
  
  arrat$Forecast<-as.integer(arrat$Forecast)


 marrat <- left_join(farrat1, xy3, by=c("Dates"))

  names(marrat)<-c("Dates","Forecast","Lower", "Upper","Cases")
  
 print (marrat)
 
 #stop(-1)
 


if(tulosta_svg==1)
{
	
	#svg(filename=tulosfilee1, width=6, height=3, pointsize=12)
	svg(filename=tulosfilee1, width=9, height=5, pointsize=12)

}

if(plottaa==1)
{



  metodi="loess"

	print ("Ggplotting ...")

#	ggplot(arrat, aes(x =Dates , y = Cases)) +ylim(yala, yyla)+xlim(as.Date(datelimits3, format="%d/%m/%Y") )+
	ggplot(marrat, aes(x =Dates , y = Forecast)) +
	
	#geom_ribbon(aes(x=Dates, ymin = Lower, ymax = Upper), fill = "lightblue") +
	#geom_line(aes(y=Forecast), size=2, color="blue")+
	
	#geom_line(aes(y=Lower), size=1)+
	#geom_line(aes(y=Upper), size=1, color="red", linetype = "dashed")+
	#annotate(geom="text", x=as.Date("2021/5/1"), y=300, size=5, label="Punainen kayra arvioitu maksimi",
    #          color="red")+
	ylim(yala, yyla) +xlim(as.Date(datelimits3, format="%d/%m/%Y") )+	
	
		
	#ggtitle("Koronavirustapauksia /pv ennuste") +
		labs(title = "Ennuste: koronavirustapauksia Suomessa",
              subtitle = "jos muutos jatkuu samaa vauhtia",
              caption = "")+
	xlab("Kuukausi") + ylab("Koronatapauksia")+
	theme(title=element_text(size=17), axis.text=element_text(size=14,face="bold"),axis.title=element_text(size=15,face="bold"))+ 
	#geom_point() +
	geom_smooth( fill="#a0a0ff",span=spanni, method=metodi, level=0.99999, size=3)+
	geom_smooth( fill="#9090ff", span=spanni,method=metodi, level=0.97) +
	geom_smooth( fill="#8a08af", span=spanni, method=metodi,level=0.8) +
	geom_point(aes(y=Cases), size=3)
	
	#geom_hline(yintercept=1.0, linetype="dashed", color = "red", size=1)
	


}


if(tulosta_svg==1)
{
dev.off()
}