BeginPackage["MonteCarlo`", "DiscreteMath`Permutations`"] Pollard::usage = "Pollard[data,rep] takes a time series of population numbers at successive times, data = {n(1), ... , n(T)}, (n(t) > 0) calculates the observed slope of log n(t+1) on log n(t), and also the slopes of rep permutations of the data using the method of Pollard et al. (1987), and returns the rank of the observed slope among all the rep+1 slopes followed by rep+1." Begin["`Private`"] Slope[data_]:=Module[{x}, Coefficient[Fit[Transpose[{Drop[data,-1],Rest[data]}],{1,x},x],x]] Pollard[data_List,rep_Integer]:=Module[{x,d,sim,dstar,xstar}, x=N[Log[data]]; d=Rest[x]-Drop[x,-1]; sim=Table[dstar=Transpose[Sort[Transpose[ {RandomPermutation[Length[d]],d}]]][[2]]; xstar=FoldList[Plus,x[[1]],dstar]; Slope[xstar],{rep}]; Flatten[{Position[Ordering[Prepend[sim,Slope[x]]],1],rep+1}]]/; VectorQ[Append[data,rep],Positive] End[] EndPackage[]