version
5.00
0
frmMain
1
1
0
0
0
7
0
4
9
0
0
0
0
0
0
0
3
0
0
0
238
268
0
0
0
Sub designer
addform(frmMain,"Normal curve","",220,220,220)@
addarraylist(frmmain,arrlstYP,10,95,80,25)@
addarraylist(frmmain,arrlstXP,10,65,80,25)@
addbutton(frmmain,btnExit,175,240,55,20,"EXIT",212,208,200,0,0,0,True,True,9)@
addradiobtn(frmmain,rbtnAbove,5,165,55,20,"Above",220,220,220,0,0,0,True,True,True,8)@
addradiobtn(frmmain,rbtnBelow,5,190,55,20,"Below",220,220,220,0,0,0,True,True,False,8)@
addradiobtn(frmmain,rbtnBetween,5,215,65,20,"Between",220,220,220,0,0,0,True,True,False,8)@
addtextbox(frmmain,tbOutsideRight,130,240,35,20,"",255,255,255,0,0,0,True,True,False,8)@
addtextbox(frmmain,tbOutsideLeft,70,240,35,20,"",255,255,255,0,0,0,True,True,False,8)@
addradiobtn(frmmain,rbtnOutside,5,240,60,20,"Outside",220,220,220,0,0,0,True,True,False,8)@
addlabel(frmmain,lblAnd,105,215,22,20,"and",220,220,220,0,0,0,True,True,8)@
addtextbox(frmmain,tbBelow,70,190,35,20,"",255,255,255,0,0,0,True,True,False,8)@
addtextbox(frmmain,tbAbove,70,165,35,20,"",255,255,255,0,0,0,True,True,False,8)@
addtextbox(frmmain,tbBetweenRight,130,215,35,20,"",255,255,255,0,0,0,True,True,False,8)@
addlabel(frmmain,lblOr,105,240,20,20,"or",220,220,220,0,0,0,True,True,8)@
addtextbox(frmmain,tbBetweenLeft,70,215,35,20,"",255,255,255,0,0,0,True,True,False,8)@
addlabel(frmmain,lblArea,50,15,140,20,"Shaded area:",220,220,220,0,0,0,True,True,8)@
addtextbox(frmmain,tbSD,195,190,35,20,"",255,255,255,0,0,0,True,True,False,8)@
addlabel(frmmain,lblSD,155,190,23,20,"SD",220,220,220,0,0,0,True,True,8)@
addtextbox(frmmain,tbMean,195,165,35,20,"",255,255,255,0,0,0,True,True,False,8)@
addlabel(frmmain,lblMean,155,165,35,20,"Mean",220,220,220,0,0,0,True,True,8)@
End Sub
@EndOfDesignText@'TO DO
'change mean and SD in graph when changed in textbox
'change as soon as textbox is entered
Sub Globals
x1 = 10
GrWidth = 215
GrHeight = 130
x2 = x1+GrWidth
y1 = 10
y2 = y1+GrHeight
Mean = 100
SD = 16
Above = Mean-1*SD
Below = Mean-1*SD
BetweenLeft = Mean-1*SD
BetweenRight = Mean+1*SD
OutsideLeft = Mean-1*SD
OutsideRight = Mean+1*SD
StepSize = 0.01
End Sub
Sub App_Start
SIP(FALSE)
Init
frmMain.Show
tbAbove.Focus = TRUE
End Sub
Sub Init
frmMain.ForeLayer = TRUE
'set rectangle for graph
frmMain.FLine(x1,y1,x2,y2,cBlack,B)
'set initial values x-axis graph
'8 marks are given on x-axis from -4*SD to +4*SD
InitXaxis
'add first point y = zeroline
arrlstXP.Add(x1+GrWidth/2-4/8*GrWidth)
arrlstYP.Add(y2)
'create grapharray in small steps
FOR c = -4 to 4 STEP StepSize
y = 1/(SQRT(2*cPI))*cE^(-0.5*c^2)
y = y*GrHeight*2
arrlstXP.Add(x1+GrWidth/2+c/8*GrWidth)
arrlstYP.Add(y2-y)
NEXT
'add last point y = zeroline
arrlstXP.Add(x1+GrWidth/2+4/8*GrWidth)
arrlstYP.Add(y2)
'draw graph
frmMain.FPolygon(arrlstXP,0,arrlstYP,0,arrlstXP.Count,cBlack)
'set initial values
tbAbove.Text = Above
tbBelow.Text = Below
tbMean.Text = Mean
tbSD.Text = SD
tbBetweenLeft.Text = BetweenLeft
tbBetweenRight.Text = BetweenRight
tbOutsideLeft.Text = OutsideLeft
tbOutsideRight.Text = OutsideRight
rbtnAbove_Click
End Sub
Sub InitXaxis
frmMain.FErase(0,y2+6,frmMain.Width,y2+16)
FOR c = -4 to 4
frmMain.FLine(x1+GrWidth/2+c/8*GrWidth,y2,x1+GrWidth/2+c/8*GrWidth,y2+5,cBlack)
frmMain.FDrawString(Mean+c*SD,6,x1+GrWidth/2+c/8*GrWidth-StrLength(Mean+c*SD)*6/2,y2+8,x1+GrWidth/2+c/8*GrWidth-StrLength(Mean+c*SD)*6/2,y2+18,cBlack)
NEXT
End Sub
Sub AreaFill(StartX,EndX,Colour)
FOR c = StartX to EndX STEP StepSize
frmMain.FLine(arrlstXP.Item((4+c)/StepSize+1),y2-1,arrlstXP.Item((4+c)/StepSize+1),arrlstYP.Item((4+c)/StepSize+1)+1,cBlue)
NEXT
End Sub
Sub AreaErase
FOR c = -3 to 3 STEP StepSize
frmMain.FErase(arrlstXP.Item((4+c)/StepSize+1),y2-1,arrlstXP.Item((4+c)/StepSize+1),arrlstYP.Item((4+c)/StepSize+1)+1)
NEXT
End Sub
Sub tbMean_GotFocus
Mean = tbMean.Text
InitXaxis
End Sub
Sub tbSD_GotFocus
SD = tbSd.Text
InitXaxis
End Sub
Sub rbtnAbove_Click
AreaErase
IF tbAbove.Text < Mean then
Value = -(tbAbove.Text-Mean)/SD
lblArea.Text = "Shaded area: "&FORMAT(NormalP(Value), "N6")
ELSE
Value = (tbAbove.Text-Mean)/SD
FORMAT(1-NormalP(Value), "N6")
END IF
AreaFill((tbAbove.Text-Mean)/SD,3,cBlue)
End Sub
Sub rbtnBelow_Click
AreaErase
IF tbBelow.Text < Mean then
Value = -(tbBelow.Text-Mean)/SD
lblArea.Text = "Shaded area: "&FORMAT(1-NormalP(Value), "N6")
ELSE
Value = (tbBelow.Text-Mean)/SD
lblArea.Text = "Shaded area: "&FORMAT(NormalP(Value), "N6")
END IF
AreaFill(-3,(tbBelow.Text-Mean)/SD,cBlue)
End Sub
Sub rbtnBetween_Click
'between = below right - below left
'positive numbers are NormalP(X)
'negative numbers have to converted to 1-NormalP(-x)
AreaErase
IF tbBetweenLeft.Text > tbBetweenRight.Text then
'switch left and right
temp = tbBetweenLeft.Text
tbBetweenLeft.Text = tbBetweenRight.Text
tbBetweenRight.Text = temp
END IF
IF (tbBetweenLeft.Text < Mean) AND (tbBetweenRight.Text < Mean) then
'left and right both < Mean
Right = -(tbBetweenRight.Text-Mean)/SD
Left = -(tbBetweenLeft.Text-Mean)/SD
lblArea.Text = "Shaded area: "&FORMAT(1-NormalP(Right)-(1-NormalP(Left)),"N6")
ELSE IF (tbBetweenLeft.Text < Mean) AND (tbBetweenRight.Text >= Mean) then
'left < Mean, right >= Mean
Right = (tbBetweenRight.Text-Mean)/SD
Left = -(tbBetweenLeft.Text-Mean)/SD
lblArea.Text = "Shaded area: "&FORMAT(NormalP(Right)-(1-NormalP(Left)),"N6")
ELSE IF (tbBetweenLeft.Text >=Mean) AND (tbBetweenRight.Text >= Mean) then
'left and right both >= Mean
Right = (tbBetweenRight.Text-Mean)/SD
Left = (tbBetweenLeft.Text-Mean)/SD
lblArea.Text = "Shaded area: "&FORMAT(NormalP(Right)-NormalP(Left),"N6")
END IF
AreaFill((tbBetweenLeft.Text-Mean)/SD,(tbBetweenRight.Text-Mean)/SD,cBlue)
End Sub
Sub rbtnOutside_Click
AreaErase
'outside = below left + above right
'below: positive numbers are NormalP(X)
'below: negative numbers are 1-NormalP(-x)
'above: positive numbers are 1-NormalP(-x)
'above: negative numbers are Normal(-x)
AreaErase
IF tbOutsideLeft.Text > tbOutsideRight.Text then
'switch left and right
temp = tbOutsideLeft.Text
tbOutsideLeft.Text = tbOutsideRight.Text
tbOutsideRight.Text = temp
END IF
IF (tbOutsideLeft.Text < Mean) AND (tbOutsideRight.Text < Mean) then
'left and right both < Mean
Left = -(tbOutsideLeft.Text-Mean)/SD
Right = -(tbOutsideRight.Text-Mean)/SD
lblArea.Text = "Shaded area: "&FORMAT(1-NormalP(Left)+NormalP(Right),"N6")
ELSE IF (tbOutsideLeft.Text < Mean) AND (tbOutsideRight.Text >= Mean) then
Left = -(tbOutsideLeft.Text-Mean)/SD
Right = (tbOutsideRight.Text-Mean)/SD
'left < Mean, right >= Mean
lblArea.Text = "Shaded area: "&FORMAT(1-NormalP(Left)+(1-NormalP(Right)),"N6")
ELSE IF (tbOutsideLeft.Text >=Mean) AND (tbOutsideRight.Text >= Mean) then
'left and right both >= Mean
Left = (tbOutsideLeft.Text-Mean)/SD
Right = (tbOutsideRight.Text-Mean)/SD
lblArea.Text = "Shaded area: "&FORMAT(NormalP(Left)+(1-NormalP(Right)),"N6")
END IF
AreaFill(-3,(tbOutsideLeft.Text-Mean)/SD,cBlue)
AreaFill((tbOutsideRight.Text-Mean)/SD,3,cBlue)
End Sub
Sub btnExit_Click
AppClose
End Sub
Sub NormalZ(X)
RETURN 1/(SQRT(2*cPI))*cE^(-0.5*X^2)
'RETURN 1/(SD*SQRT(2*cPI))*cE^(-0.5*((X-Mean)/SD)^2)
End Sub
Sub NormalP(A)
' Returns P(A) for the Standard Normal Distribution as defined by
' Abramowitz & Stegun. This is the Probability that a value is less
' than A, i.e. Area under the curve defined by NormalZ to the left of A
' Only handles values A >= 0 otherwise exception raised.
' Accuracy: Absolute Error < 7.5e-8 }
B1 = 0.319381530
B2 = -0.356563782
B3 = 1.781477937
B4 = -1.821255978
B5 = 1.330274429
T = 1 / (1 + 0.2316419 * A)
T2 = T^2
T3 = T2*T
T4 = T2^2
T5 = T4*T
RETURN (1.0 - NormalZ(A) * (B1*T + B2*T2 + B3*T3 + B4*T4 + B5*T5))
End Sub
'If the amount people spend per day in a given store is Normally Distributed with a mean of $200 and a Standard Deviation of $25, then what is the probability that a given person will spend at most $250?
'
'Now our above Approximation doesn't reference the Mean or the Standard Deviation. This is because the above all refers to the Standard Normal Distribution which has a Mean of 0 and a Standard Deviation of 1.
'
'We can easily transform any other Normal Distribution Problem into a Standard Normal Distribution problem by using:
'
'
'
'This however can result in Negative values, and our above Approximation only works for Non-negatives. We can get around this by using the Symmetry of the graph, resulting in the following routine:
'
'function NormalDistP (const Mean, StdDev, AVal: Extended): Single;
'{Returns the Probability of (X < AVal) for a Normal Distribution
' with given Mean and Standard Deviation.
' Standard Deviation must be > 0 or function will result in an
' exception.
' Accuracy: Absolute Error < 7.5e-8 }
'var
' Z: Extended;
' Lower: Boolean;
'begin
' if FloatIsZero (StdDev) or (StdDev < 0) then
' raise EMathError.Create ('Standard Deviation must be positive')
' else
' begin
' Z := (AVal - Mean) / StdDev; // Convert to Standard (z) value
' Lower := Z < 0; // If Negative use Symmetry to calculate
' if Lower then
' Z := (Mean - AVal) / StdDev;
' Result := NormalP (Z); // Access function
' if Lower then // If Negative use Symmetry to calculate
' Result := 1.0 - Result;
' end;
'end;
'
'
'With the above we could solve the specified problem by simply calling:
'
' Value := NormalDistP (200, 25, 250);
'
'But what do we do if we want a Probability that the spend at least $300 or perhaps we want the probability that they spend between $150 and $250?
'
'Once again we use the Symmetry of the above Curve and some simple Algebra.
'
'
'
'function NormalDistQ (const Mean, StdDev, AVal: Extended): Single;
'{ Returns the Probability of (X > AVal) for a Normal Distribution
' with given Mean and Standard Deviation.
' Standard Deviation must be > 0 or function will result in an
' exception.
' Accuracy: Absolute Error < 7.5e-8 }
'begin
' Result := 1 - NormalDistP (Mean, StdDev, Aval);
'end;
'
'
'Thus the at least $300 problem would be done by simply calling:
'
' Value := NormalDistQ (200, 25, 300);
'
'
'
'function NormalDistA (const Mean, StdDev, AVal, BVal: Extended): Single;
'{ Returns the Probability of (AVal < X < BVal) for a Normal Distribution
' with given Mean and Standard Deviation.
' Standard Deviation must be > 0 or function will result in an
' exception.
' Accuracy: Absolute Error < 7.5e-8 }
'begin
' if SameFloat (AVal, BVal) then
' Result := 0
' else
' begin
' Result := NormalDistP (Mean, StdDev, BVal) -
' NormalDistP (Mean, StdDev, AVal);
' if BVal < AVal then
' Result := -1 * Result;
' end;
'end;
'
'
'Thus the at between $150 and $250 problem would be done by simply calling:
'
' Value := NormalDistA (200, 25, 150, 250);
'
'and notice we have constructed it so that the following would give the same result:
'
' Value := NormalDistA (200, 25, 250, 150);
'Normal Distributions in Delphi - Part II
'This issue we continue our look at developing Statistical Routines to use in your Delphi applications, and in particular conclude looking at using Normal distributions. These will be designed to use open array parameters where possible so that you can use them for standard arrays or with the dynamic arrays of Delphi 4 and later. While our freeware ESBMaths unit contains many stats routines (and more being added), if you want a good commercial stats library for Delphi, I would recommend either our ESBPCS or Turbo Power's SysTools 3.
'
'In a latter issue, we will look at Math and Stats resources available on the Internet.
'
'Reference
'Formulae and other information used in this article rely on Handbook of Mathematical Functions with Formulas, Graphs and Mathematical Tables by Milton Abramowitz and Irene A. Stegun, Dover. ISBN 0-486-61272-4. While not for the faint hearted this is a worthwhile book for Programmers interested in Mathematics.
'
'Another Problem
'After much research, it is found that the average life of the GloBrite Light Bulb under continuous use is 75 days, with a standard deviation of 15 days. GloBrite offers a free replacement if a Light Bulb expires within 50 days of purchase. What percentage of Light Bulbs would expire resulting in a free replacement?
'
'Given the routines we encountered in Part 1 we could easily write code to handle the above:
'
' Value := NormalDistP (75, 15, 50);
'
'Which would tell us that about 4.78% of the light bulbs sold would get a free replacement.
'
'However, what if decreased production costs now mean that GloBrite is prepared to replace 7.5% of all bulbs sold under their replacement policy - to what could they extend their period for free replacement?
'
'In other words, we want to be able to do the reverse of the above.
'
'Inverse Normal Distribution
'When it comes to probability distributions we normally need routines to calculate the Probabilities as we did last issue and we also need routines to calculate the Inverses. Most Approximations for Inverses of Distributions are only accurate to a few decimal places, though this is often more than sufficient.
'
'As mentioned last issue, Abramowitz & Stegun provides many different Approximations that can be used and in particular the following is one for the Inverse Normal Distribution adapted into Delphi:
'
'function InvNormalQ (const P: Extended): Extended;
'const
' C0: Extended = 2.515517;
' C1: Extended = 0.802853;
' C2: Extended = 0.010328;
' D1: Extended = 1.432788;
' D2: Extended = 0.189269;
' D3: Extended = 0.001308;
'var
' T: Extended;
'begin
' T := Sqrt (Ln (1.0 / Sqr (P)));
' Result := T - (C0 + C1 * T + C2 * Sqr (T)) /
' (1.0 + D1 * T + D2 * Sqr (T) + D3 * Sqr (T) * T);
'end;
'
'
'The above works for probabilities that are greater than 0 and less than or equal to 0.5. The resultant error is less than 0.00045 - thus giving us 3 decimal place accuracy.
'
'We now add a routine with a more useful interface and Exception Handling:
'
'{ Returns the Original AVal for Pr (X < AVal) = PVal if Less is
' True, or Pr (X > AVal) = PVal if Less is False,
' for a Normal Distribution with given Mean and Standard Deviation.
'
' Standard Deviation must be > 0 or function will return an Exception.
' PVal (Probability) must be 0 < PVal < 1 or function will return an Exception.
'
' Approximation has an absolute error < 4.5 x 10^4
' Thus 3 figure accuracy.
'
'}
'function InvNormalDist (const Mean, StdDev, PVal: Extended; const Less: Boolean): Extended;
'var
' P: Extended;
' Lower: Boolean;
'begin
' if not FloatIsPositive (StdDev) then
' raise EMathError.Create ('Standard Deviation must be Positive');
'
' if (PVal <= 0) or (PVal >= 1) then
' raise EMathError.Create ('Probability must be strictly between 0 and 1');
'
' if Less then
' begin
' P := 1 - PVal;
' Lower := P > 0.5;
' if Lower then
' P := PVal;
' end
' else
' begin
' P := PVal;
' Lower := P > 0.5;
' if Lower then
' P := 1 - PVal;
' end;
' Result := InvNormalQ (P);
'
' if Lower then
' Result := Mean - (Result * StdDev)
' else
' Result := Mean + (Result * StdDev);
'end;
'
'
'Thus we can now answer the problem that was proposed:
'
' Value := InvNormalDist (75, 15, 0.75, True);
'
'Which results in finding out that the replacement period could be extended to 53 days.
'
'
'MathRtns and Demo
'From this Issue we now have a Zip file called MathsCorner.Zip - that will be updated each issue. It contains a Unit called MathRtns.Pas which contains all the routines so far covered in the Maths Corner as well as some extras.
'
'You will also find NormDist2 Project which demonstrates the above problem.