Blackouts.q
From InterSciWiki
blackouts <- c(7500000,2100000,2085000,1660000,1600000,1500000,1500000,1500000,1300000,1140000,899000,877000,875000,870000,725000,660000,650000,600000,600000,598000,570000,557354,500000,500000,500000,490000,464000,460000,404000,400000,400000,385000,375000,363476,360000,350000,350000,320831,315000,312000,300000,300000,300000,300000,290000,284000,272000,272000,258000,257718,250000,246000,246000,242910,240000,238000,235000,234000,230000,219000,210882,207200,206000,203000,200000,191000,190000,173000,173000,166000,164500,163000,160000,160000,160000,158000,148000,147000,146000,145000,145000,142000,133000,130000,130000,128000,126000,124000,122000,120000,120000,115000,114500,114000,113200,112000,106850,100000,100000,100000,100000,95630,95000,94285,92000,92000,91000,90000,88000,82500,81000,80000,80000,75000,75000,74000,71000,71000,71000,70000,70000,70000,70000,66005,65000,65000,63500,62000,60000,60000,60000,60000,60000,59000,58000,56000,56000,55000,55000,53000,51000,51000,50462,50000,50000,50000,
+ 50000,50000,50000,50000,50000,50000,48000,46000,45000,43696,43000,40911,40000,40000,40000,39500,38500,37000,36073,35000,33000,32000,32000,30500,30001,30000,29900,29000,29000,26334,25000,25000,25000,25000,25000,25000,24506,24000,20000,20000,19000,18819,18351,18000,18000,17000,15000,15000,14273,12000,11529,11000,10300,10000,10000,9000,8000,7500,5300,4150,2900,2000,1800,1646,1000) blackouts.min <- 230000 blackouts.tailprob <- sum(blackouts>=blackouts.min)/length(blackouts) blackouts.pareto <- pareto.fit(blackouts,blackouts.min) blackouts.random <- blackouts.min*((1-ppareto(blackouts,blackouts.min,blackouts.pareto$exponent,lower.tail=FALSE))^(1/(1-blackouts.pareto$exponent))) blackouts.tsal <- tsal.fit(blackouts,xmin=1000) #blackouts.tsal <- tsal.fit(ifelse(blackouts<=blackouts.min,blackouts,blackouts.random),xmin=1000) plot.survival.loglog(blackouts,ylab="Cumulative probability",main="Survival function of blackouts data", sub="zagged=empirical, solid=Tsallis, dashed=Pareto, ") curve(ptsal(x,blackouts.tsal$shape,blackouts.tsal$scale,lower.tail=FALSE),add=TRUE,col="blue") curve(blackouts.tailprob*ppareto(x,blackouts.min,blackouts.pareto$exponent,lower.tail=FALSE),from=blackouts.min,col="red",lty="dashed",add=TRUE) blackouts.tsal $type [1] "tsal" $q [1] 1.565541 $kappa [1] 113940.9 $shape [1] 1.768219 $scale [1] 201472.4 $loglike [1] -2788.139 $n [1] 211 $xmin [1] 1000 $method [1] "mle.equation"
