1974年1月から2014年8月までの太陽黒点数


SSN:monthly averages , DEV:standard deviation

Sys.time()
## [1] "2014-09-23 18:55:04 JST"

Package

#install.packages("knitr",dependencies=T)
#install.packages("tseries",dependencies=T)
#install.packages("quantmod",dependencies=T)
#install.packages("ggplot2",dependencies=T)
library(knitr)
library(tseries)
library(quantmod) # http://cran.r-project.org/web/packages/quantmod/quantmod.pdf
library(ggplot2)

Read Data

url<-"http://solarscience.msfc.nasa.gov/greenwch/spot_num.txt"
ss<-read.table(url,header=T,skip=0,stringsAsFactor=FALSE,as.is=T)
ss$date<-as.Date(paste(ss$YEAR,"/",ss$MON,"/1",sep=""))
ss<-ss[,c(5,3,4)]
str(ss)
## 'data.frame':    3188 obs. of  3 variables:
##  $ date: Date, format: "1749-01-01" "1749-02-01" ...
##  $ SSN : num  58 62.6 70 55.7 85 83.5 94.8 66.3 75.9 75.5 ...
##  $ DEV : num  24.1 25.1 26.6 23.6 29.4 29.2 31.1 25.9 27.7 27.7 ...
summary(ss)
##       date                 SSN             DEV      
##  Min.   :1749-01-01   Min.   :  0.0   Min.   : 0.0  
##  1st Qu.:1815-05-24   1st Qu.: 15.7   1st Qu.:11.3  
##  Median :1881-10-16   Median : 42.2   Median :19.2  
##  Mean   :1881-10-15   Mean   : 52.1   Mean   :20.2  
##  3rd Qu.:1948-03-08   3rd Qu.: 76.5   3rd Qu.:26.8  
##  Max.   :2014-08-01   Max.   :253.8   Max.   :90.2
head(ss,3)
##         date  SSN  DEV
## 1 1749-01-01 58.0 24.1
## 2 1749-02-01 62.6 25.1
## 3 1749-03-01 70.0 26.6
tail(ss,3)
##            date  SSN  DEV
## 3186 2014-06-01 71.0 28.7
## 3187 2014-07-01 72.5 45.0
## 3188 2014-08-01 74.7 20.8
adf.test(ss$SSN)
## 
##  Augmented Dickey-Fuller Test
## 
## data:  ss$SSN
## Dickey-Fuller = -7.036, Lag order = 14, p-value = 0.01
## alternative hypothesis: stationary
summary(diff(ss$SSN))
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  -75.10   -8.30    0.00    0.01    8.40   94.90

Plot

qplot(SSN,DEV,data=ss)

plot of chunk unnamed-chunk-4

qplot(date,SSN,data=ss,geom=c("line","point"))

plot of chunk unnamed-chunk-4

qplot(ss$date[-1],diff(ss$SSN))

plot of chunk unnamed-chunk-4

ggplot(ss,aes(x=SSN))+geom_histogram(binwidth=10,fill="grey",alpha=0.5,color="black")

plot of chunk unnamed-chunk-4

acf(ss$SSN,type="correlation",main="ACF of Sunspot Number",lag=240,ci=c(0.9,0.95))

plot of chunk unnamed-chunk-4