S-PLUS
The statistical methods described here use S-Plus, a statistical package available at http://www.mathsoft.com/splus.html.
Additional software and extensions for S-Plus are available at http://lib.stat.cmu.edu/S/.
The command S-Plus will start up the package on most platforms. This will automatically create a .Data directory where the S-Plus data will be stored. The examples shown in these pages assume that whatever files are being scanned in are in the same directory where S-Plus was started, as are any object files which are necessary. If this is not the case, the path should be included in the scan calls. The C files which are used in some of the methods above need to be compiled to an object file, which is what the S-Plus routines actually use.
Splus is a statistical package which allows the user a wide flexibility in usage. Simple Splus objects include not only numbers, vectors and matrices, but also inherent and user-defined functions. Following is a sample session, with comments denoted by "#". ">" is the prompt. Splus # Starts Splus from shell X11() #Starts graphics window. #All of the usual arithmetic and logical operators are present. 5+3 > 8 #The "help()" function allows for on-line help about any Splus object, be it #function, etc. help(X11) #Provides help on the command X11(). #To attach a directory different form .Data, use the attach command: attach("/usr/joe/.Data") #For more information help(attach) #One can create and dimension a vector of length n, filled with zeros, by #the command "vect <- rep(0,n)". "<-" or "_" is the assignment operator. hestimates <- rep(0,50) #This is a nice way to create a vector which can be used later, e.g. in a #FOR loop. The vector elements are referred to as "hestimates[i]". #One can simulate stable series with the command "rstab". sample.series <- rstab(1000, index=1.5) # Length 1000, alpha=1.5. tsplot(sample.series) # Produces time series plot. #One can also simulate non-symmetric stable series -- the parametrization #used in this case is special -- see Samorodnitsky and Taqqu, Chapman #and Hall, 1994, p. 43. sample.series[1] > 0.8194 lsfit(sample.series)$coefficients > 1 .5 # $coefficients produces a subset of the output, namely, the coefficients # of a least square fit(Intercept, slope). # (Residuals, etc. can also be obtained). stem(sample.series) # Produces stem plot (median, quartiles, etc.) # One can simulate long-range dependent series (FARIMA) with the command # "arima.fracdiff.sim" . For the right syntax, use help: help(arima.fracdiff.sim) #or look at the function: arima.fracdiff.sim sample.series2 <- arima.fracdiff.sim(model=list(d=0.3),n=1000) sample.series2[1] > 1.293 #Most of the statistical methods take a vector as the "data" variable, and #run the vector through some manipulations, and usually come up with #numerical results and a plot. For example, #(1) Create and dimension an array: farima.data_array(0,dim=c(50,1000)) # Creates an array of zeros 50 by 1000. Members of array can be referenced as # farima.data[i,j]. A whole column or row of array can be referenced as # farima.data[i,] or farima.data[,j]. #(2) Fill in array # Now we put simulated long-range dependent series into array, making 50 # series, each length 1000. We use a FOR loop. for(i in 1:50){ farima.data[i,]<-arima.fracdiff.sim(model=list(d=0.1),n=1000)} #(3) Use an estimation procedure: # Suppose we want to estimate the scaling parameter H using the Aggregated # Variance method. We will use the function "plotvar" contained in the # file "file.S". The source command reads in the contents of this file. source("file.S") #Now we estimate H using our Aggregated Variance method. for (i in 1:50}{ hestimates[i]_plotvar(farima.data[i,],plot=False)} #(4) Examine results: #Since we are analyzing 50 series, we do not usually want to see plots. #Now one can find statistics for the estimator of H. mean(hestimates) sqrt(var(hestimates)) hist(hestimates) #returns histogram of this vector. #(5) Writing results to outside file, sink all output: sink("file_name") mean(hestimates) sqrt(var(hestimates)) sink() # closes file. # To see what objects one has defined, use the command: objects() # This lists everything which has been defined and not deleted. # To delete things use rm(). rm(sample.series) # To see what an object is, type in the name: plotvar #returns the definition of the function plotvar hestimates #returns the values of the vector hestimates #To split graphics window into 4 sub-plots, 2 by 2. par(mfrow=c(2,2)) #For info on axes, title, captions, etc. help(plot) #To save a graphic (plot): One has to regenerate the output and direct it #to a file instead of the X11 window. Last device created and not yet #closed is the default. Open file: postscript("name_of_file.ps") plotvar(sample.series, plot=T) #generate plot (portrait is default). dev.off() #close file #Alternatively, the print button on the X11() window saves a postscript #file as a by-product of printing. (Saved in directory where S-Plus was #started. The default is landscape, but portrait printing can be done by #pressing the "print orientation" button in the X11() window. # Control-C interrupts. Control-Z suspends. #Quit q()