Skip to content

Instantly share code, notes, and snippets.

@amerberg amerberg/weibull-sum.R
Last active Aug 29, 2015

What would you like to do?
#set seed for replicability
#generate 10^3 deviations, 10^6 replicates
w = replicate(10^6, rweibull(10^3, 0.9, 1))
#recall Weibull is subexponential when shape parameter is < 1.
#animate histogram of ratios of max to sum for first n
sums = rep(0, 10^6)
maxs = rep(0, 10^6)
increment = 10
for (n in increment*(1:(1000/increment))) {
start = n-increment+1
end = n
#recalculate sum/max for each trial by adding `increment` new terms
sums = sums + colSums(w[start:end,])
maxs = colMaxs(rbind(w[start:end,], maxs))
print(qplot(maxs/sums, geom="histogram", main=sprintf("Histogram of Max/Sum (n=%d)",n),
xlab="ratio", xlim=c(0,1), ylim=c(0,3e5), binwidth=0.001))
}},"weibull.gif", interval=0.2, ani.width=1000, ani.height=1000)
#Image shows ratio concentrated about 0 as number of terms increases to 1000.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.