Code: Select all
!npop = 1000 'size of the population
wfcreate u !npop 'create an unstructured workfile
!mu = 5 'population mean
!sig = 0.5 'population standard deviation
series ypop = !mu + !sig*@rnorm 'hypothetical population
!nrep = 100 'number of draws
!nsamp = 50 'sample size
!alpha = 0.05 'significance level
%insmpl = @pagesmpl 'population range
%outsmpl = "1 "+@str(!nsamp) 'sample range
'Series objects are bound by the frequency of the workfile
if !nrep>!npop then
@uiprompt("Number of draws should be less than the population size!")
!nrep = 100
endif
smpl 1 !nrep 'adjust the page sample
series meanval 'create a series to hold the mean of each sample
series lbound 'create a series to hold the lower bounds of the confidence intervals
series ubound 'create a series to hold the upper bounds of the confidence intervals
group confidence 'create a group to store the confidence intervals
confidence.add ubound lbound meanval 'add series to the group object
%false="" 'create a blank string
'Loop for drawing samples
for !i=1 to !nrep
ypop.resample(outsmpl=%outsmpl,insmpl=%insmpl,permute,suffix=ysamp) 'sampling
meanval(!i) = @mean(ysamp) 'mean values
lbound(!i) = @mean(ysamp) - @qnorm(1-!alpha/2)*@stdev(ysamp)/@sqrt(!nsamp) 'lower bound
ubound(!i) = @mean(ysamp) + @qnorm(1-!alpha/2)*@stdev(ysamp)/@sqrt(!nsamp) 'upper bound
if !mu>ubound(!i) or !mu<lbound(!i,2) then
%false=%false+" "+@str(!i) 'false positives (type I errors)
endif
statusline !i of !nrep 'monitor the progress
next
%confchart = @getnextname("confchart") 'name the chart
!rotate=0 '=1 to rotate the chart
if !rotate=1 then
freeze({%confchart}) confidence.errbar(rotate) 'generate the chart
{%confchart}.draw(line,bottom,color(red),width(1)) 5 'draw the population mean
for !j=1 to @wcount(%false)
!s = @val(@word(%false,!j))
{%confchart}.draw(shade,left,rgb(230,184,183)) !s 'highlight false positives
next
else
freeze({%confchart}) confidence.errbar
{%confchart}.draw(line,left,color(red),width(1)) 5
for !j=1 to @wcount(%false)
!s = @val(@word(%false,!j))
{%confchart}.draw(shade,bottom,rgb(230,184,183)) !s
next
endif
%header = "# samples = "+@str(!nrep)+" | # rejects = "+@str(@wcount(%false))+" | sample size = "+@str(!nsamp)+" | moments of sampling distribution = ("+@str(@mean(meanval),"f.2")+","+@str(@stdevp(meanval),"f.3")+")"
{%confchart}.addtext(t) {%header}
{%confchart}.setupdate(m) 'enables the slider bar
show {%confchart} 'display the chart