Question: Trouble Finding Standard Deviation for Density Function

I am trying to find the standard deviation for a piecewise density function (f).  The interval is (0, sqrt(3))  I am convinced the density function is what I want.  I calculate the expected value of the density function (evE) and the answer is correct -- about 0.66145

When I try to calculate the standard deviation (stdE), I get an answer that is "off" by a large degree.  Via simulating values, I should get a value of about 0.24936.  Each time I "re-run" the calculation, I get varying results, all of which are "off" by a large degree.

I am only guessing, but my integration function might be missing some sort of assumption and/or option.

My work is attached.  Does anyone know what I am doing wrong?

restart

with(Statistics)

f := -t*piecewise(t <= 0, 0, t <= 1, t*(t^3+6*Pi*t-8*t^2-4*Pi), t <= sqrt(2), -(2*(4*t^2*sqrt(t^2-1)*arcsin((t^2-2)/t^2)+4*t^2*arctan(sqrt(t^2-1))*sqrt(t^2-1)+(t^4+(2*Pi+3)*t^2-4*t*Pi+3*Pi-1/2)*sqrt(t^2-1)-8*t^4+4*t^2+4))/sqrt(t^2-1), t <= sqrt(3), (8*arctan((t^3+t^2-3*t-1)/(sqrt(t^2-2)*(t^3-t^2-t-1)))*t*sqrt(t^2-2)-4*arccot((t^2-t-1)/sqrt(t^2-2))*t^2*sqrt(t^2-2)+4*t*sqrt(t^2-2)*(t+2)*arctan((t^2+t-1)/sqrt(t^2-2))+(8*t^2-4)*sqrt(t^2-2)*arcsin((t^2-3)/(t^2-1))+16*arcsin(sqrt(t^2-2)/sqrt(t^2-1))*sqrt(t^2-2)-16*arcsin(1/sqrt(t^2-1))*sqrt(t^2-2)+((t^2+5)*sqrt(t^2-2)-8*t^2+16)*(t^2+1))/sqrt(t^2-2), sqrt(3) < t, 0)

-t*piecewise(t <= 0, 0, t <= 1, t*(t^3+6*Pi*t-8*t^2-4*Pi), t <= 2^(1/2), -2*(4*t^2*(t^2-1)^(1/2)*arcsin((t^2-2)/t^2)+4*t^2*arctan((t^2-1)^(1/2))*(t^2-1)^(1/2)+(t^4+(2*Pi+3)*t^2-4*t*Pi+3*Pi-1/2)*(t^2-1)^(1/2)-8*t^4+4*t^2+4)/(t^2-1)^(1/2), t <= 3^(1/2), (8*arctan((t^3+t^2-3*t-1)/((t^2-2)^(1/2)*(t^3-t^2-t-1)))*t*(t^2-2)^(1/2)-4*arccot((t^2-t-1)/(t^2-2)^(1/2))*t^2*(t^2-2)^(1/2)+4*t*(t^2-2)^(1/2)*(t+2)*arctan((t^2+t-1)/(t^2-2)^(1/2))+(8*t^2-4)*(t^2-2)^(1/2)*arcsin((t^2-3)/(t^2-1))+16*arcsin((t^2-2)^(1/2)/(t^2-1)^(1/2))*(t^2-2)^(1/2)-16*arcsin(1/(t^2-1)^(1/2))*(t^2-2)^(1/2)+((t^2+5)*(t^2-2)^(1/2)-8*t^2+16)*(t^2+1))/(t^2-2)^(1/2), 3^(1/2) < t, 0)

(1)

plot(f, t = 0 .. sqrt(3))

 

evE := simplify(int(t*f, t = 0 .. sqrt(3)))

(1/105)*(294*(2^(1/2)-3/2)*(-2+3^(1/2))*ln(1+2^(1/2))+168*(2^(1/2)-3/2)*(-2+3^(1/2))*ln(1+3^(1/2))-252*(2^(1/2)-3/2)*(-2+3^(1/2))*arccoth(2^(1/2))+((-14*Pi-84*ln(2)-19)*3^(1/2)+28*Pi+168*ln(2)+50)*2^(1/2)+(21*Pi+126*ln(2)+20)*3^(1/2)-42*Pi-252*ln(2)-58)/((2*2^(1/2)-3)*(-2+3^(1/2)))

(2)

evalf(evE)

.6617071838

(3)

stdE := simplify(sqrt(int(f*(t-evE)^2, t = 0 .. sqrt(3))))

(1/210)*(-4808607312140100*(2^(1/2)-3880899/2744210)*(3^(1/2)-70226/40545)*ln(1+2^(1/2))^2+(-5495551213874400*(2^(1/2)-3880899/2744210)*(3^(1/2)-70226/40545)*ln(1+3^(1/2))+8243326820811600*(2^(1/2)-3880899/2744210)*(3^(1/2)-70226/40545)*arccoth(2^(1/2))+((457962601156200*Pi+2747775606937200*ln(2)+631290043332900)*3^(1/2)-793214493249360*Pi-4759286959496160*ln(2)-1093426429714584)*2^(1/2)+(-647656921614780*Pi-3885941529688680*ln(2)-892778941072368)*3^(1/2)+1121774694224184*Pi+6730648165345104*ln(2)+1546338486359352)*ln(1+2^(1/2))-1570157489678400*(2^(1/2)-3880899/2744210)*(3^(1/2)-70226/40545)*ln(1+3^(1/2))^2+(4710472469035200*(2^(1/2)-3880899/2744210)*(3^(1/2)-70226/40545)*arccoth(2^(1/2))+((261692914946400*Pi+1570157489678400*ln(2)+360737167618800)*3^(1/2)-453265424713920*Pi-2719592548283520*ln(2)-624815102694048)*2^(1/2)+(-370089669494160*Pi-2220538016964960*ln(2)-510159394898496)*3^(1/2)+641014110985248*Pi+3846084665911488*ln(2)+883621992205344)*ln(1+3^(1/2))-3532854351776400*(2^(1/2)-3880899/2744210)*(3^(1/2)-70226/40545)*arccoth(2^(1/2))^2+(((-392539372419600*Pi-2355236234517600*ln(2)-541105751428200)*3^(1/2)+679898137070880*Pi+4079388822425280*ln(2)+937222654041072)*2^(1/2)+(555134504241240*Pi+3330807025447440*ln(2)+765239092347744)*3^(1/2)-961521166477872*Pi-5769126998867232*ln(2)-1325432988308016)*arccoth(2^(1/2))+((-392539372419600*ln(2)^2+(-130846457473200*Pi-180368583809400)*ln(2)-10903871456100*Pi^2-30061430634900*Pi+1205966071986042)*3^(1/2)+679898137070880*ln(2)^2+(226632712356960*Pi+312407551347024)*ln(2)+18886059363080*Pi^2+52067925224504*Pi-2088794509063308)*2^(1/2)+(555134504241240*ln(2)^2+(185044834747080*Pi+255079697449248)*ln(2)+15420402895590*Pi^2+42513282908208*Pi-1705493574764511)*3^(1/2)-961521166477872*ln(2)^2+(-320507055492624*Pi-441810996102672)*ln(2)-26708921291052*Pi^2-73635166017112*Pi+2954001523727854)^(1/2)*(5*2^(1/2)-7)*(-5+3*3^(1/2))/((70*2^(1/2)-99)*(15*3^(1/2)-26)*(-2+3^(1/2))*(4*3^(1/2)-7)*(12*2^(1/2)-17)*(2*2^(1/2)-3))

(4)

evalf(stdE)

7031348.295

(5)
 

NULL

Download StdE.mw

Please Wait...