Carl Love

Carl Love

28035 Reputation

25 Badges

12 years, 319 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are replies submitted by Carl Love

@Wes Good question. It helps to visualize the p-value as the area in the "tails" of the distribution. That might be the left tail, the right tail, or the sum of both tail areas; it depends on what type of inequality is in your "alternative hypothesis"--- x < a (left), x > a (right), or x <> a (both). Here, x is some statistical parameter (such as mean) that you're testing, and a is some fixed real number. (Let me know if you need more details on any of that.)

Suppose that T is the test distribution (such as Chi^2 or Student's t) (which is not the underlying distribution of the raw data, such as Normal), and t is the value of the test statistic (based on some computation, often elaborate, done with the raw data). Then the p-value for a 

  • right-tailed test is 1 - CDF(T, t),
  • left-tailed test is CDF(T, t),
  • two-tailed test on a symmetric test distribution (such as Student's t) is
    2*CDF(T, -abs(t)),
  • two-tailed test on an asymmetric distribution (such as Chi^2 or F)... I need to defer to another expert's opinion on that; I don't know off the top of my head.

@nm When a symbol appears as a left operand of ->, that is its declaration.

@Carl Love I wrote above:

  • It's fairly easy to customize the number of terms in the series to your value of Digits. The series is alternating, so Leibniz's Theorem on the error of truncated alternating series applies.

Actually, it's so easy to do this that I couldn't resist writing it up and posting it here. Using a procedure of only two lines to numerically evaluate the asymptotic series, we can get results that are several times more efficient and several digits more accurate than the default numeric evaluation of your function. 

In the code below, I've omitted your factor of 4
 

restart:

f:= x-> x*ln(1+1/x);

proc (x) options operator, arrow; x*ln(1+1/x) end proc

The next line of code gives a formal presentation of the asymptotic series of f(x). This
is not part of the computation. It's only here to serve as a guide to coding the procedure.

simplify(eval(convert(f(1/x), FormalPowerSeries), x= 1/x));

Sum((-x)^(-n)/(n+1), n = 0 .. infinity)

An alternative way of viewing the series, also just a guide:

value(eval(%, infinity= 6));

1-(1/2)/x+(1/3)/x^2-(1/4)/x^3+(1/5)/x^4-(1/6)/x^5+(1/7)/x^6

This is the procedure for the numeric evaluation of f(x). There are two main factors to its efficiency:

1. 

No explicit powers of -x need to be computed. Instead, the power is computed as the previous term's
power divided by -x.

2. 

There's no need to explicitly strive for a particular accuracy. Instead, we just keep adding terms until
they no longer make a difference.


This procedure requires 1D Input.

`evalf/f1`:= (x::float)->  local r, r1:= 1, p:= 1, n:= 1;
    if x>1 then do r:= r1 until (r1+= (p/= -x)/++n) = r else f(x) fi
:

f1:= x-> `if`(x::float, evalf('f1'(x)), f(x)):

Digits:= 16:
seq(f1(10.^k), k= 1..20);

.9531017980432486, .9950330853168082, .9995003330835331, .9999500033330833, .9999950000333330, .9999995000003333, .9999999500000033, .9999999950000000, .9999999995000000, .9999999999500000, .9999999999950000, .9999999999995000, .9999999999999500, .9999999999999950, .9999999999999995, 1, 1, 1, 1, 1

That procedure is many times more efficient than the default numeric evaluation:

X:= evalf(['rand'() $ 2^14]):

forget(evalf); gc(); CodeTools:-Usage(f1~(X)):

memory used=21.25MiB, alloc change=0 bytes, cpu time=109.00ms, real time=107.00ms, gc time=0ns

forget(evalf); gc(); CodeTools:-Usage(f~(X)):

memory used=72.20MiB, alloc change=0 bytes, cpu time=469.00ms, real time=400.00ms, gc time=109.38ms

And it's also many times more accurate:

forget(evalf);

evalf[32](evalf[16]([f,f1](99.)) - evalf[32]([f,f1](99.)));

[-0.997717133689829526e-14, 0.2282866310170378e-16]

 

Download AsymptEvalf.mw

@Carl Love Here's an easy-to-understand and purely algebraic alternative:

subs(x= x-1, collect(subs(x= x+1, p), x));

The first-degree term expands automatically. To prevent that, you can do

subs(x= ``(x-1), collect(subs(x= x+1, p), x));

 

@MathPrincess123 Your concept of long division is completely wrong. The first term of the quotient should be the quotient of the highest-degree terms, not the lowest-degree terms. Perhaps you are being confused (visually) by your sort order. Descending degrees is the usual order, just like ordinary integers (where the degrees are the exponents of 10).

Also, it can be clearly seen that my program doesn't make any "change" to the polynomials other than replacing q with Q, because there are no fractional or negative exponents. The quotient is the same as it would be under any division process, whether by hand or computer. In this particular case, the quotient of the highest-degree terms is the degree-0 term of the quotient.

This does not need to be your "last question". Feel free to keep asking until you completely understand.

@MathPrincess123 

I thought that you might want to handle the negative-exponent case, but I didn't include it because the quotient and remainder need to be treated separately. To restore the original variable, the exponents of the remainder need to be reduced in the negative direction, but the exponents of the quotient do not.

So, add this procedure:

QuoRem:= proc(n, d, x::name)
local X, ND, R, m, q, r;
    (ND,R):= FractionalExponentAdjuster([n,d], x, X);
    m:= min(0, ldegree~(ND, X)[]);
    q:= quo(expand(ND/~X^m)[], X, r);
    #Remainder gets downshifted, but quotient does not:
    R(q), R(expand(X^m*r))
end proc
:

To use it, simply do

(qNSR, rNSR):= QuoRem(NSR, deltaa, q);
expand(qNSR*deltaa + rNSR - NSR);
#Test---should get 0

And, no, it was not "an easy question".

@Wes I don't think that you understand the meaning of CDFPDF, and Quantile. If is a continuous random variable, and x is a real number, then CDF(X, x) is the probability that X <= x. So, CDF returns a probability, a number between 0 and 1. PDF(X, x) is the derivative with respect to x of CDF(X, x). It is mainly of theoretical interest; there's no practical value to evaluating it at a specific number such as with your attempt PDF(Y, 0.05)Quantile(X, p) is the functional inverse of CDF(X, x); in other words, it's the solution x of the equation CDF(X, x) = p. Note that p here is a probability.

@max125 What output were you expecting? If it's [2, 7], it can be obtained as 

(f@op)~([[2,0], [3,4]]); 

@vv After doing several experiments replacing the bindings with modules (either nested or unnested), I agree with you that it's a bug; the result should be (a+b)*a. For example, nested modules:

restart:
A:= module() export a:= :-a+b, B:= module() export b:= a - :-b; end module; end module:
use A in use B in a*b end use end use;

 

@vv My guess is that the order that the automatic simplifications are applied makes the difference.

@vv It agrees, like this:

`%%%%subs`(
    {x= 0, y= 0}, 
    `%%%subs`(a= a+b+x, `%%subs`(%subs(a= a+b+x, b= a-b+y), a*b))
):
(value@@4)(%);
                                   2
                            (a + b) 

(I realize that there are problems with this example!)

@vv I think that you're ignoring that the outer use also changes the a in b= a-b.

I make sense of this example via this simulation:

`%%%subs`(a= a+b, `%%subs`(%subs(a= a+b, b= a-b), a*b));
  %%%subs(a = a + b, %%subs(%subs(a = a + b, b = a - b), a b))

value(%);
              %%subs(a = a + b, %subs(b = a, a b))

value(%);
                           /            2\
                      %subs\a = a + b, a /

value(%);
                                   2
                            (a + b) 

 

@The function Show a worksheet for what you did get. B.A is definitely a 3x3 matrix, both in Maple and in math.

@The function The angle bracket constructor for row vectors uses the vertical bar character (|) as the element separator:

A:= <1 | 3 | 4>;

You can also enter it as a column vector and take its transpose:

A:= <1, 3, 4>^%T;

B:= <5, -2, 1>;

Multiplication of vectors and matrices is non-commutative (i.e., order matters), so it uses a different operator: the dot `.` rather than `*`.

A.B; B.A;

In this case, A.B will be a single number, and B.A will be a 3x3 matrix.

@ecterrab You wrote:

  • ...unwith unloads the package and forget clears remember tables from results obtained with the package, both concepts apply to any package. This is expected and what one needs to do, @Carl Love....

Look again at my 2nd example, repeated below. The result that I wanted to remove hadn't been obtained with the package, and, had I been a naive user, I would have had no reason to suspect that the result was different simply because Physics had been loaded (yet hadn't been otherwise used).

#In this example, no command from Physics is explicitly used, yet the results
#from top-level :-diff are different.
restart:	
with(Physics): 	
:-diff(:-conjugate(f(x)), x);

Of course, in actual practice there's a lengthy computation that involves diff rather than a simple call to diff as shown. Thus, a naive user would have no way of knowing what needs to be forgotten.

First 73 74 75 76 77 78 79 Last Page 75 of 708