Nothing
#################################################################################
##
## R package pmoments by Alexios Ghalanos Copyright (C) 2008
## This file is part of the R package pmoments.
##
## The R package pmoments is free software: you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation, either version 3 of the License, or
## (at your option) any later version.
##
## The R package pmoments is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
## GNU General Public License for more details.
##
#################################################################################
#----------------------------------------------------------------------------------
# Partial Moments Optimization plots
# Weights Plot
weightsPlot = function(object, col = NULL, legend = TRUE, title=NULL)
{
UseMethod("weightsPlot")
}
.weightsPlot=function(object, col = NULL, legend = TRUE, title=NULL)
{
if (is.null(col)) col = rainbow(ncol(object@weights))
weights = object@weights
pos.weights = +0.5 * (abs(weights) + weights)
neg.weights = -0.5 * (abs(weights) - weights)
ymax = max(rowSums(pos.weights))
ymin = min(rowSums(neg.weights))
range = ymax - ymin
ymax = ymax + 0.005 * range
ymin = ymin - 0.005 * range
dim = dim(weights)
range = dim[1]
xmin = 0
xmax = range + 0.2 * range
if(!legend){
barplot(t(pos.weights), space = 0, ylab = "",
ylim = c(ymin, ymax), col = col, border = "grey")
} else {
legendtext = colnames(object@weights)
if(is.null(legendtext)){
for(i in 1:dim[2]){legendtext[i] = paste("Asset", i, sep = " ")}
}
barplot(t(pos.weights), space = 0, ylab = "", xlim = c(xmin, xmax),
ylim = c(ymin, ymax), col = col, border = "grey")
legend("topright", legend = legendtext, bty = "n", cex = 0.8,
fill = col)
}
barplot(t(neg.weights), space = 0, add = TRUE, col = col, border = "grey")
targetRisk = object@riskMeasure
targetReturn = object@rewardMeasure
targetRaR = targetReturn/targetRisk
nSigma = length(targetRisk)
nLabels = 6
M = c(0, ( 1:(nSigma %/% nLabels) ) ) *nLabels + 1
nPrecision = 3
axis(1, at = M, labels = signif(targetRisk[M], nPrecision),cex.axis=0.8)
axis(3, at = M, labels = signif(targetReturn[M], nPrecision),cex.axis=0.8)
mtext("Target Risk", side = 1, line = 2, cex = 0.9)
mtext("Target Return", side = 3, line = 2, cex = 0.9)
mtext("Weights (sum)", side = 2, line = 2, cex = 0.9)
lines(x = c(0, nSigma), c(1, 1), col = "grey", lty = 3)
lines(x = c(0, nSigma), c(0, 0), col = "grey", lty = 3)
minIndex = which.min(targetRisk)
minRisk = signif(min(targetRisk), 3)
# Add vertical Line at max risk/return:
tR = which(targetRaR==max(targetRaR))
tRR = signif(targetRisk[tR], 3)
abline(v = tR, col = "steelblue", lty = 2, lwd = 2)
mtext(paste(
object@solverType, "|", object@solver, "|", "minRisk =", minRisk),
side = 4, adj = 0, col = "grey", cex = 0.7)
mtext(paste(
object@solverType, "|", object@solver, "|", "optimalRisk =", tRR),
side = 4, adj = 0, padj=1.2 , col = "grey", cex = 0.7)
if(is.null(title)) mtext("Weights", adj = 0, line = 2.5, font = 2, cex = 0.8) else mtext(title, adj = 0, line = 2.5, font = 2, cex = 0.8)
box()
invisible()
}
setMethod("weightsPlot", signature(object="PMSOLVER"), .weightsPlot)
weightsPie = function(object, legend = TRUE, title=NULL)
{
UseMethod("weightsPie")
}
.weightsPie<-function(object, legend = TRUE, title=NULL)
{
if(length(object@riskMeasure)!=1) stop("function not available for frontier")
weights<-(object@weights)
wsort<-sort.int(weights,index.return = TRUE)
weights=weights[1,wsort$ix]
labelsx = names(weights)
labels = names(weights)
Radius = 1
if(length(weights)>15){
z=which(abs(weights)<0.01)
labels[z]=""
Radius = 0.85
}
negw<-which(weights<0)
posw<-which(weights>=0)
negc<-.colorgradient(n = length(negw), low.col = 0.0, high.col=0.1, saturation = 1)
posc<-.colorgradient(n = length(posw), low.col = 0.6, high.col=0.65, saturation = 1)
mycol=vector(mode="character",length=length(weights))
mycol=c(negc,posc)
par(no.readonly = TRUE)
pie(abs(weights), labels = labels, col = mycol, radius = Radius,cex=0.6)
box()
# Add Title:
if(is.null(title)){
title(main = paste("Partial Moment Optimization Weights [moment: ", object@moment,"]",sep=""))
} else{
title(main=title)
}
# Add Info:
mtext(paste("Solver Type :", object@solverType, " | reward measure :", round(object@rewardMeasure,5)," | risk measure :",round(object@riskMeasure,5),sep=""),
side = 4, adj = 0, col = "grey", cex = 0.7)
mtext(paste("Blue Tones: Positive Weights | Red Tones: Negative Weights",sep=""),
side = 1, adj=0, col = "grey", cex = 0.7)
mtext(paste("Assets with abs(weight)<1% do not display labels",sep=""),
side = 1, adj=0, padj=1.5, col = "grey", cex = 0.7)
mtext(paste("pmoments",sep=""),
side = 2, adj=0,padj=-1, col = "grey", cex = 0.7)
if(legend){
legendWeights = as.character(round(100*weights, digits = 1))
legendWeights = paste(legendWeights, "%",sep=" ")
legend(1.1,1, legend = labelsx, bty = "n", cex = 0.8,
fill = mycol)
legend(1.3,1, legend = legendWeights,bty = "n", cex = 0.8)
}
# Return Value:
invisible()
}
setMethod("weightsPie", signature(object="PMSOLVER"), .weightsPie)
# Frontier Plots & Lines
.plot.pmSolver<-function(x,...)
{
if(length(x@riskMeasure)==1) stop("x is not a solver type IV (frontier) object")
risk<-x@riskMeasure
reward<-x@rewardMeasure
riskFree=x@riskFree
tangent = x@tangent
xlabel<-paste("Lower Partial Moment [m = ",x@moment,"]",sep="")
ylabel<-"Expected Return"
minrisk = which(risk==min(risk))
plot(risk[1:(minrisk-1)],reward[1:(minrisk-1)],type="p",col="black",ylab=ylabel,xlab=xlabel,xlim=c(0.5*min(risk),1.1*max(risk)),
ylim=c(0,1.05*max(reward)),lwd=2,...)
points(risk[minrisk:length(risk)],reward[minrisk:length(reward)],col="grey",lwd=2)
points(risk[tangent],reward[tangent],col="steelblue",lwd=3)
points(risk[minrisk],reward[minrisk],col="violet",lwd=3)
title("LPM Frontier Plot",col="black",lty=2)
leg.txt=c("Frontier(>minRisk)", "Frontier(<minrisk)",paste("Tangency(riskFree=",round(riskFree,4),")",sep=""), "minRisk")
legend(x="bottomright", yjust=0,legend=leg.txt,lwd=3, lty=c(1,1,1,3), col=c('black','grey','steelblue','violet'), cex=0.8)
mtext(paste("pmoments package"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.7)
grid()
invisible()
}
.lines.pmSolver<-function(x,...)
{
if(length(x@riskMeasure)==1) stop("x is not a solver type IV (frontier) object")
risk<-x@riskMeasure
reward<-x@rewardMeasure
lines(risk,reward,...)
invisible()
}
.points.pmSolver<-function(x,...)
{
if(length(x@riskMeasure)==1) stop("x is not a solver type III (frontier) object")
risk<-x@riskMeasure
reward<-x@rewardMeasure
points(risk,reward,...)
invisible()
}
setMethod("plot", signature(x="PMSOLVER",y="missing"), .plot.pmSolver)
setMethod("lines", signature(x="PMSOLVER"), .lines.pmSolver)
setMethod("points", signature(x="PMSOLVER"), .points.pmSolver)
#----------------------------------------------------------------------------------
# Partial Moments Utility Plots
.plot.pmu<-function(x,...)
{
dx<-sort(x@data)
dy<-x@utility
plot(dx,dy,type="l",ylab="Utility",xlab="X",...)
leg.text<-paste("LM=",x@lmoment,", UM=",x@umoment,", LC=",x@lrc,", UC=",x@urc,sep="")
#legend("topleft",legend=leg.text,text.col=,fill=1,cex=0.7,inset=c(0.03,0.06))
title("Upper-to-Lower Partial Moment Utility Function")
grid()
invisible()
}
.lines.pmu<-function(x,...)
{
op <- options()
options(warn = -1)
dx<-sort(x@data)
dy<-x@utility
lines(dx,dy,...)
options(op)
invisible()
}
setMethod("plot", signature(x="PMU",y="missing"), .plot.pmu)
setMethod("lines", signature(x="PMU"), .lines.pmu)
#----------------------------------------------------------------------------------
# Partial Moments Expectation Plots
# Sub-methods for single/multiple data type and estimation method
.plot.pme<-function(x,...)
{
if(x@datatype=="single"){
mm=x@method
for(i in 1:length(mm)){
sol<-switch(mm[i],
boot=.sbootplot(x,...),
jack=.sjackplot(x,...),
sample=.ssampleplot(x,...),
analytical=.sanalyticalplot(x,...))
}} else{
munq=rev(unique(x@method))
for(i in 1:length(munq)){
sol<-switch(munq[i],
boot=.mbootplot(x,...),
jack=.mjackplot(x,...),
sample=.msampleplot(x,...),
analytical=.manalyticalplot(x,...))
}}
invisible()
}
setMethod("plot", signature(x="PME",y="missing"), .plot.pme)
#----------------------------------------------------------------------------------
# Single DataType Plots
.sbootplot<-function(x,...)
{
n = length(x@moment)
cl=rainbow(n=n)
if(n>16){
scr=floor(n/16)
z=n-16*floor(n/16)
start=dev.next(which = dev.cur())
tally=0
for(j in 1:scr){
dev.new(start+j)
par(mfrow=c(4,4))
for(i in 1:16){
tmp=x@pm$boot[[i+tally]]$thetastar
hist(tmp,col=cl[i+tally],ylab="frequency",xlab=paste("moment = ",round(x@moment[i+tally],2),sep=""),main="")
mtext(paste("pmoments : bootstrap"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("Boostrapped Partial Moment Histograms [",x@datanames,"]","\n(page...",j,")",sep=""), outer=TRUE, line=-2, cex=0.75)
tally=tally+16
}
if(z!=0){
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(4,4))
for(i in 1:z){
tmp=x@pm$boot[[i+tally]]$thetastar
hist(tmp,col=cl[i+tally],ylab="frequency",xlab=paste("moment = ",round(x@moment[i+tally],2),sep=""),main="")
mtext(paste("pmoments : bootstrap"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("Boostrapped Partial Moment Histograms [",x@datanames,"]","\n(page...",j,")",sep=""), outer=TRUE, line=-2, cex=0.75)
}
} else{
d=.divisortable(n)
par(mfrow=c(d[1],d[2]))
for(i in 1:n){
tmp=x@pm$boot[[i]]$thetastar
hist(tmp,col=cl[i],ylab="frequency",xlab=paste("moment = ",round(x@moment[i],2),sep=""),main="")
mtext(paste("pmoments : bootstrap"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.5)
mtext(paste("mean: ",round(mean(tmp),3*(as.integer(x@moment[i]/2)))," | sd: ",round(sd(tmp),3*(as.integer(x@moment[i]/2))),sep=""),
side = 4, adj = 0, padj=1.2, col = "gray", cex = 0.5)
title(paste("Boostrapped Partial Moment Histograms [",x@datanames,"]",sep=""),outer=TRUE,cex=0.85, line=-2)
}
}
invisible()
}
.sjackplot<-function(x,...)
{
n = length(x@moment)
cl=rainbow(n=n)
if(n>16){
scr=floor(n/16)
z=n-16*floor(n/16)
start=dev.next(which = dev.cur())
tally=0
for(j in 1:scr){
dev.new(start+j)
par(mfrow=c(4,4))
for(i in 1:16){
tmp=x@pm$jack[[i+tally]]$jack.values
hist(tmp,col=cl[i+tally],ylab="frequency",xlab=paste("moment = ",round(x@moment[i+tally],2),sep=""),main="")
mtext(paste("pmoments : JackKnife"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("JackKnife Partial Moment Histograms [",x@datanames,"]","\n(page...",j,")",sep=""), outer=TRUE, line=-2, cex=0.75)
tally=tally+16
}
if(z!=0){
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(4,4))
for(i in 1:z){
tmp=x@pm$jack[[i+tally]]$jack.values
hist(tmp,col=cl[i+tally],ylab="frequency",xlab=paste("moment = ",round(x@moment[i+tally],2),sep=""),main="")
mtext(paste("pmoments : JackKnife"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("JackKnife Partial Moment Histograms [",x@datanames,"]","\n(page...",j,")",sep=""), outer=TRUE, line=-2, cex=0.75)
}
} else{
d=.divisortable(n)
par(mfrow=c(d[1],d[2]))
for(i in 1:n){
tmp=x@pm$jack[[i]]$jack.values
hist(tmp,col=cl[i],ylab="frequency",xlab=paste("moment = ",round(x@moment[i],2),sep=""),main="")
mtext(paste("pmoments : JackKnife"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.6)
mtext(paste("JackKnife s.e.(",round(x@pm$jack[[i]]$jack.se,5),")",sep=""), side = 4, adj = 0, padj=1.2, col = "gray", cex = 0.6)
title(paste("JackKnife Partial Moment Histograms [",x@datanames,"]",sep=""),outer=TRUE,cex=0.85, line=-2)
}
}
invisible()
}
.ssampleplot<-function(x,...)
{
n = length(x@moment)
cl=rainbow(n=n)
if(n>16){
scr=floor(n/16)
z=n-16*floor(n/16)
start=dev.next(which = dev.cur())
tally=0
tmpdata=x@data
threshold=x@threshold
tmpsdata=sort(tmpdata)
for(j in 1:scr){
dev.new(start+j)
par(mfrow=c(4,4))
for(i in 1:16){
tmp=.pmValues(tmpsdata, threshold=threshold, moment=x@moment[i+tally], type=x@tail, standardize=x@standardize)
if(x@tail=="lower"){
pd=max(which(tmpsdata<=threshold))
zt=1:pd
plot(tmpsdata[zt],tmp[zt],col=cl[i+tally],ylab=paste("moment = ",round(x@moment[i+tally],2)), xlab="data",main="",cex.lab=0.9)
} else{
pd=min(which(tmpsdata>threshold))
zt=pd:length(tmpsdata)
plot(tmpsdata[zt],tmp[rev(zt)],col=cl[i+tally],ylab=paste("moment = ",round(x@moment[i+tally],2)), xlab="data",main="",cex.lab=0.9)
}
plot(tmpsdata[zt],tmp[zt],col=cl[i+tally],ylab=paste("moment = ",round(x@moment[i+tally],2)), xlab="data",main="",cex.lab=0.9)
mtext(paste("pmoments : sample"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("Cumulative Sample Partial Moment-to-Threshold [",x@datanames,"]","\n Tail=",x@tail," (page...",j,")",sep=""), outer=TRUE, line=-2, cex=0.75)
tally=tally+16
}
if(z!=0){
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(4,4))
for(i in 1:z){
tmp=.pmValues(tmpsdata, threshold=threshold, moment=x@moment[i+tally], type=x@tail, standardize=x@standardize)
if(x@tail=="lower"){
pd=max(which(tmpsdata<=threshold))
zt=1:pd
plot(tmpsdata[zt],tmp[zt],col=cl[i+tally],ylab=paste("moment = ",round(x@moment[i+tally],2)), xlab="data",main="",cex.lab=0.9)
} else{
pd=min(which(tmpsdata>threshold))
zt=pd:length(tmpsdata)
plot(tmpsdata[zt],tmp[rev(zt)],col=cl[i+tally],ylab=paste("moment = ",round(x@moment[i+tally],2)), xlab="data",main="",cex.lab=0.9)
}
mtext(paste("pmoments : sample"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("Cumulative Sample Partial Moment-to-Threshold [",x@datanames,"]","\n Tail=",x@tail," (page...",j,")",sep=""), outer=TRUE, line=-2, cex=0.75)
}
} else{
d=.divisortable(n)
tmpdata=x@data
threshold=x@threshold
tmpsdata=sort(tmpdata)
par(mfrow=c(d[1],d[2]))
for(i in 1:n){
tmp=.pmValues(tmpsdata, threshold=threshold, moment=x@moment[i], type=x@tail, standardize=x@standardize)
if(x@tail=="lower"){
pd=max(which(tmpsdata<=threshold))
zt=1:pd
plot(tmpsdata[zt],tmp[zt],col=cl[i],ylab=paste("moment = ",round(x@moment[i],2)), xlab="data",main="",cex.lab=0.9)
} else{
pd=min(which(tmpsdata>threshold))
zt=pd:length(tmpsdata)
plot(tmpsdata[zt],tmp[rev(zt)],col=cl[i],ylab=paste("moment = ",round(x@moment[i],2)), xlab="data",main="",cex.lab=0.9)
}
mtext(paste("pmoments : sample"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
title(paste("Cumulative Sample Partial Moment-to-Threshold [",x@datanames,"]\n Tail=",x@tail,sep=""),outer=TRUE,cex=0.85, line=-2)
}
}
invisible()
}
.sanalyticalplot<-function(x,...)
{
# sample versus analytical
# 2 plots (density and sample with partial moment shading on threshold
# and barplot with overlayed moments (sample vs analytical)
print("single dataset/analytical method: no plot method for this type of data")
}
#----------------------------------------------------------------------------------
.mbootplot<-function(x,...)
{
# we need to distinguish the methods:
zz=which(x@method=="boot")
n = length(zz)
datanames=x@datanames[zz]
cl=rainbow(n=n)
if(n>16){
scr=floor(n/16)
z=n-16*floor(n/16)
#start=dev.new(dev.next(which = dev.cur())+1)
tally=0
for(j in 1:scr){
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(4,4))
for(i in 1:16){
tmp=x@pm[[zz[i+tally]]]$thetastar
hist(tmp,col=cl[zz[i+tally]],ylab="frequency",xlab=paste("moment = ",round(x@moment[zz[i+tally]],2),sep=""),main=datanames[zz[i+tally]])
mtext(paste("pmoments : bootstrap"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("Boostrapped Partial Moment Histograms (page...",j,")",sep=""), outer=TRUE, line=-1, cex=0.75)
tally=tally+16
}
if(z!=0){
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(4,4))
for(i in 1:z){
tmp=x@pm[[zz[i+tally]]]$thetastar
hist(tmp,col=cl[zz[i+tally]],ylab="frequency",xlab=paste("moment = ",round(x@moment[zz[i+tally]],2),sep=""),main=datanames[zz[i+tally]])
mtext(paste("pmoments : bootstrap"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("Boostrapped Partial Moment Histograms (page...",scr+1,")",sep=""), outer=TRUE, line=-1, cex=0.75)
}
} else{
d=.divisortable(n)
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(d[1],d[2]))
for(i in 1:n){
tmp=x@pm[[zz[i]]]$thetastar
hist(tmp,col=cl[i],ylab="frequency",xlab=paste("moment = ",round(x@moment[zz[i]],2),sep=""),main=datanames[i])
mtext(paste("pmoments : bootstrap"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.5)
mtext(paste("mean: ",round(mean(tmp),3*(as.integer(x@moment[zz[i]]/2)))," | sd: ",round(sd(tmp),3*(as.integer(x@moment[zz[i]]/2))),sep=""),
side = 4, adj = 0, padj=1.2, col = "gray", cex = 0.5)
title("Boostrapped Partial Moment Histograms", outer=TRUE,cex=0.85, line=-1)
}
}
}
.mjackplot<-function(x,...)
{
# we need to distibuish the methods:
zz=which(x@method=="jack")
n = length(zz)
datanames=x@datanames[zz]
cl=rainbow(n=n)
if(n>16){
scr=floor(n/16)
z=n-16*floor(n/16)
#start=dev.next(which = dev.cur())
tally=0
for(j in 1:scr){
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(4,4))
for(i in 1:16){
tmp=x@pm[[zz[i+tally]]]$jack.values
hist(tmp,col=cl[zz[i+tally]],ylab="frequency",xlab=paste("moment = ",round(x@moment[zz[i+tally]],2),sep=""),main=datanames[zz[i+tally]])
mtext(paste("pmoments : JackKnife"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("JackKnife Partial Moment Histograms (page...",j,")",sep=""), outer=TRUE, line=-1, cex=0.75)
tally=tally+16
}
if(z!=0){
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(4,4))
for(i in 1:z){
tmp=x@pm[[zz[i+tally]]]$jack.values
hist(tmp,col=cl[zz[i+tally]],ylab="frequency",xlab=paste("moment = ",round(x@moment[zz[i+tally]],2),sep=""),main=datanames[zz[i+tally]])
mtext(paste("pmoments : JackKnife"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("JackKnife Partial Moment Histograms (page...",scr+1,")",sep=""), outer=TRUE, line=-1, cex=0.75)
}
} else{
d=.divisortable(n)
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(d[1],d[2]))
for(i in 1:n){
tmp=x@pm[[zz[i]]]$jack.values
hist(tmp,col=cl[zz[i]],ylab="frequency",xlab=paste("moment = ",round(x@moment[zz[i]],2),sep=""),main=datanames[zz[i+tally]])
mtext(paste("pmoments : JackKnife"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.5)
mtext(paste("mean: ",round(mean(tmp),3*(as.integer(x@moment[zz[i]]/2)))," | sd: ",round(sd(tmp),3*(as.integer(x@moment[zz[i]]/2))),sep=""),
side = 4, adj = 0, padj=0.8, col = "gray", cex = 0.5)
mtext(paste("JackKnife s.e.(",round(x@pm[[zz[i]]]$jack.se,5),")",sep=""), side = 4, adj = 0, padj=1.2, col = "gray", cex = 0.6)
title("JackKnife Partial Moment Histograms", outer=TRUE,cex=0.85, line=-2)
}
}
}
#----------------------------------------------------------------------------------
# used in choosing the mfrow
.divisortable<-function(n)
{
z=matrix(c(1,1,1,
2,2,1,
3,2,2,
4,2,2,
5,2,3,
6,2,3,
7,2,4,
8,2,4,
9,3,3,
10,3,4,
11,3,4,
12,3,4,
13,4,4,
14,4,4,
15,4,4,
16,4,4),ncol=3,byrow=TRUE)
d=which(n==z[,1])
return(z[d,2:3])
}
#----------------------------------------------------------------------------------
.msampleplot<-function(x,...)
{
# we need to distinguish the methods:
zz=which(x@method=="sample")
n = length(zz)
datanames=x@datanames[zz]
cl=rainbow(n=n)
if(n>16){
scr=floor(n/16)
z=n-16*floor(n/16)
#start=dev.new(dev.next(which = dev.cur())+1)
tally=0
for(j in 1:scr){
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(4,4))
for(i in 1:16){
tmpdata=x@data[,zz[i+tally]]
threshold=x@threshold[zz[i+tally]]
tmpsdata=sort(tmpdata)
tmp=.pmValues(tmpsdata, threshold=threshold, moment=x@moment[zz[i+tally]], type=x@tail[1], standardize=x@standardize[zz[i+tally]])
if(x@tail[1]=="lower"){
pd=max(which(tmpsdata<=threshold))
zt=1:pd
plot(tmpsdata[zt],tmp[zt],col=cl[zz[i+tally]],ylab=paste("moment = ",round(x@moment[zz[i+tally]],2)), xlab="data",main=datanames[i+tally],cex.lab=0.9)
} else{
pd=min(which(tmpsdata>threshold))
zt=pd:length(tmpsdata)
plot(tmpsdata[zt],tmp[rev(zt)],col=cl[zz[i+tally]],ylab=paste("moment = ",round(x@moment[zz[i+tally]],2)), xlab="data",main=datanames[i+tally],cex.lab=0.9)
}
mtext(paste("pmoments : bootstrap"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("Cumulative Sample Partial Moment-to-Threshold (page...",j,")",sep=""), outer=TRUE, line=-1, cex=0.75)
tally=tally+16
}
if(z!=0){
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(4,4))
for(i in 1:z){
tmpdata=x@data[,zz[i+tally]]
threshold=x@threshold[zz[i+tally]]
tmpsdata=sort(tmpdata)
tmp=.pmValues(tmpsdata, threshold=threshold, moment=x@moment[zz[i+tally]], type=x@tail[1], standardize=x@standardize[zz[i+tally]])
if(x@tail[1]=="lower"){
pd=max(which(tmpsdata<=threshold))
zt=1:pd
plot(tmpsdata[zt],tmp[zt],col=cl[zz[i+tally]],ylab=paste("moment = ",round(x@moment[zz[i+tally]],2)), xlab="data",main=datanames[i+tally],cex.lab=0.9)
} else{
pd=min(which(tmpsdata>threshold))
zt=pd:length(tmpsdata)
plot(tmpsdata[zt],tmp[rev(zt)],col=cl[zz[i+tally]],ylab=paste("moment = ",round(x@moment[zz[i+tally]],2)), xlab="data",main=datanames[i+tally],cex.lab=0.9)
}
mtext(paste("pmoments : bootstrap"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("Cumulative Sample Partial Moment-to-Threshold (page...",scr+1,")",sep=""), outer=TRUE, line=-1, cex=0.75)
}
} else{
d=.divisortable(n)
dev.new(dev.next(which = dev.cur())+1)
par(mfrow=c(d[1],d[2]))
for(i in 1:n){
tmpdata=x@data[,i]
threshold=x@threshold[i]
tmpsdata=sort(tmpdata)
tmp=.pmValues(tmpsdata, threshold=threshold, moment=x@moment[i], type=x@tail[1], standardize=x@standardize[i])
if(x@tail[1]=="lower"){
pd=max(which(tmpsdata<=threshold))
zt=1:pd
plot(tmpsdata[zt],tmp[zt],col=cl[i],ylab=paste("moment = ",round(x@moment[i],2)), xlab="data",main=datanames[i],cex.lab=0.9)
} else{
pd=min(which(tmpsdata>threshold))
zt=pd:length(tmpsdata)
plot(tmpsdata[zt],tmp[rev(zt)],col=cl[i],ylab=paste("moment = ",round(x@moment[i],2)), xlab="data",main=datanames[i],cex.lab=0.9)
}
mtext(paste("pmoments : bootstrap"), side = 4, adj = 0, padj=0, col = "gray", cex = 0.4)
}
title(paste("Cumulative Sample Partial Moment-to-Threshold",sep=""), outer=TRUE, line=-1, cex=0.75)
}
}
.manalyticalplot<-function(x,...)
{
# sample versus analytical
# 2 plots (density and sample with partial moment shading on threshold
# and barplot with overlayed moments (sample vs analytical)
print("multiple dataset/analytical method: no plot method for this type of data")
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.