Carl Love

Carl Love

28025 Reputation

25 Badges

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

MaplePrimes Activity


These are answers submitted by Carl Love

The Maple equivalent of the Matlab syntax Q(:, 1), where Q is a Matrix, is Q[.., 1].

To answer your overall question---to wit:

  • I would like to ask if there is any command in Maple while solving the equations that are responsible for possible approximations when solving equations

If you use the command solve or any command in its "family", which include the package SolveTools and numerous other commands and packages, and your input does not include any decimal numbers, then there are no approximations used. If the command cannot find an exact solution, then it returns nothing.

What makes you think that the solution that you got is wrong?

In the worksheet attached to your original Question, you give the following formula, which you attribute to VV:

  • The total number of set partitions of n distinct items into e sets of length L, with n = e*L, is n!/e!/(L!)^e.

To be more precise, that is the number of ways of doing it if the subsets are unlabeled. If they are labeled (as in red team, green team, orange team), then the formula is n!/L!^e. Both of these formulas can be easily generalized to the situation where the subsets aren't necessarily equally sized.

Let be a list of positive integers representing the sizes of the blocks of the partition. If the blocks are labeled, then the number of partitions is add(B)!/mul(B!~). This can also be computed in Maple as combinat:-multinomial(add(B), B[]). If the blocks are unlabeled, then the number of partitions is add(B)!/mul(B!~)/mul(rhs~(Statistics:-Tally(B))!~). This can also be computed in Maple as Number(Iterator:-SetPartitionFixedSize(B)).

For B = [5,5,5], the labeled-block formula gives 756756 and the unlabeled 126126. The former leads to the probability 1/1001 and the latter 6/1001. The denominator from my original Answer is binomial(15,5)*binomial(10,5)/3!, which also equals 126126.

 

If you're expecting the result 2, as if differentiating with respect to x(t), do like this:

f:= 2*x(t)+5;
Physics:-diff(f, x(t));

I never use 2D input, so I can't help with that aspect of doing it.

It's not a bug (although it's very often majorly inconvenient), and it's not anything specific to IntegrationTools:-Change. Changing assumptions makes new copies of variables which are not automatically substituted into the expressions containing the old copies. The way to get that substitution to happen (temporarily) is to use assuming instead of the second assume and assuming additionally instead of additionally, like this:

restart;
assume(a>0,b>a,c>0,t>0);
F := Int(sqrt(d-a*c^2*t),d=0...infinity);
IntegrationTools[Change](F,-a*t*c^2+d=-y,y) assuming b > a;

or

restart;
assume(a>0,b>a,c>0,t>0);
F := Int(sqrt(d-a*c^2*t),d=0...infinity);
IntegrationTools[Change](F,-a*t*c^2+d=-y,y) assuming additionally, b > a;

 

As you probably know, there is a one-to-one correspondence between partitions of n of length k and partitions of n whose maximum element is k. Using that principle, here's a procedure that's much more efficient than my prior Answer:

PartitionFixedSize:= proc(n::posint, k::depends(integer[1..n]))
local p, P;
   [seq(1 +~ `+`(seq([0$(k-p), 1$p], p= P)), P= combinat:-partition(n-k, min(n-k, k)))]
end proc:

 

When you want to use an index on a function's name as if it were a function parameter, the index must be extracted inside the procedure using op(procname). Unfortunately, this makes the syntax a bit awkward, but there's no need to use a for loop for this. Thus, this entire section of your code:

g[1](r);
ux[1] := (x, r) -> epsilon[1][1]*x + g[1](r);
ur[1] := r -> A[1]*r + B[1]*1/r;
varepsilon[1][1] := epsilon11c;
varepsilon[1][2] := r -> (A[1]*r + B[1]*1/r)*1/r;
varepsilon[1][3] := r -> diff(ur[1](r), r);
varepsilon[1][3](R[2]);

for i from 2 to N - 1 do 
g[i](r); 
ux[i] := (x, r) -> epsilon[i][1]*x + g[i](r); 
ur[i] := r -> A[i]*r + B[i]*1/r; 
varepsilon[i][1] := epsilon11c; 
varepsilon[i][2] := r -> (A[i]*r + B[i]*1/r)*1/r; 
varepsilon[i][3] := r -> diff(ur[i](r), r); 
varepsilon[i][2](r); i;
end do;

should be replaced by

ux:= (x,r)-> epsilon[(local i:= op(procname))][1]*x + g[i](r):
ur:= r-> A[(local i:= op(procname))]*r + B[i]/r:
varepsilon:= proc(r) 
local k:= op(procname), i:= op([0,1], procname);
   if k=1 then epsilon11c
   elif k=2 then ur[i](r)/r
   elif k=3 then diff(ur[i](r), r)
   fi
end proc:   

And you would invoke it, for example, like this:

varepsilon[2][1](); varepsilon[2][2](r); varepsilon[2][3](r);

So, note that that's a slight change in the calling syntax for varepsilon[i][1](), because the parentheses are needed.

Here's one way to do it:

PartitionFixedSize:= (n::posint, k::depends(integer[1..n]))-> 
  subsindets(Statistics:-Tally~(combinat:-composition(n,k)), `=`, `$`@op)
:
PartitionFixedSize(10,2);
            {[1, 9], [2, 8], [3, 7], [4, 6], [5, 5]}

This is not the most computationally efficient way to do it, because I start with the compositions and discard the duplicates. Of course, it'd be more efficient to not generate duplicates in the first place. But my procedure has the benefit of brevity.

There's a command that's supposed to do it efficiently, but it seems completely broken. It's called Iterator:-PartitionFizedSize. The command 

{seq([seq(p)], p= Iterator:-PartitionFixedSize(10,2))}

should return the same thing, but all that it actually returns is {[9,1]}. Hopefully @Joe Riel will fix this bug soon.

In October 2018, I Posted an article and code for Fisher's exact test. This is a good situation to use it. (The other three tests mentioned in that Post don't apply to this situation, and, unlike Fisher's, they're only approximations anyway.) The contingency table to test (i.e., the observed table) is 

SurvivorTeams:= <
   5, 0, 4; 
   0, 5, 1
>:

where each column represents one of the three new teams and each row represents one of the two original teams. (It makes no difference to the final result whether the roles of the rows and columns are switched, i.e., if we were to transpose the matrix, although my code is a bit more computationally efficient if you use the more-columns-than-rows form.)

The null hypothesis is that the 15 players were assigned to the three teams independently of which of the two original teams they were from, i.e., randomly. The alternative hypothesis is that it wasn't done randomly. Like all (two-tailed[*1]) statistical hypothesis tests, this one does not and cannot give a definitive yes-or-no answer to that. Rather, they all return a p-value: the probability that the observed data (such as our matrix above)---or any other arrangement of data at least as deviant or extreme---would be observed if the null hypothesis were true. (I cannot overemphasize the importance of the underlined phrase in the previous sentence. Make sure that you understand it. It's fundamental to all statistical hypothesis testing. In a discrete case, such as the current one, deviant or extreme could be reworded as improbable.) Thus, a small value of provides evidence that the alternative hypothesis is true.

Running my procedure on the given matrix:

infolevel[ContingencyTableTests]:= 1:
interface(prettyprint= 1):
ContingencyTableTests(SurvivorTeams, method= Fisher): 
% = evalf(%);
25 tables evaluated, of which 6 are at least as improbable as the observed input.
                             6                   
                            ---- = 0.005994005994
                            1001                 

Note that this is numerically the same probability as I gave in my off-the-top-of-my-head prior Answer, although the computational process that my procedure uses is much more elaborate and automatic (it having computed the probabilities for 25 different arrangements and summing 6 of those).

Furthermore, I object to Kitonum's Reply being referred to (by mmcdara) as a "correction" of my prior Answer.

In the code below[*2], I give 5 weblink references. Of these, I most highly recommend the Mathworld article "Fisher's Exact Test"[*3], which contains this formula (totally reworded by me): The probabality of any particular matrix T 's occurence in is (&*!(R)*&*!(C)/n!)/&*!(T) where &*! is the product of the factorials of the entries (of a vector or matrix), is the row-sum vector, is the column-sum vector, n is the total of all entries, and S is the sample space of all matrices of nonnegative integers that have the same R and C. For the table SurvivorTeams, that's

(6!*9!*5!^3/15!)/(5!*5!*4!)
     1 / 1001

Or, more programmatically, that's

`&*!`:= mul@factorial~:
Pr:= proc(T::Matrix(nonnegint))
uses AT=ArrayTools;
   &*!(AT:-AddAlongDimension(T,1))*&*!(AT:-AddAlongDimension(T,2))/add(T)!/&*!(T)
end proc:

Pr(SurvivorTeams);

There are 25 possible tables with the same R and C, with 6 of those corresponding to the 6 possible permutations of the columns or 6 possible assignment of colors to the new teams, each with equal probability, and none of the other 19 have lower probability. Hence, again we arrive at the answer 6 / 1001. (The interested reader should try to construct a few of the 19 tables omitted from the sum, as Fisher did long before the computer era.)
 

[*1]two-tailed statistical hypothesis test is one where the alternative hypothesis is the direct negation of the null hypothesis. It is called a two-sided alternative. A one-sided alternative involves a directional inequality.

[*2]Here's my code, which has been substantially updated since my October 2018 Post:

#This code uses features new to Maple 2019.

ContingencyTableTests:= proc(
	O::Matrix(nonnegint), #contingency table of observed counts 
	{method::identical(Pearson, Yates, G, Fisher, Fisher[count]):= 'Pearson'}
)
description 
	"Returns p-value for Pearson's (w/ or w/o Yates's continuity correction)" 
	" or G chi-squared or Fisher's exact test."
	" All computations are done in exact arithmetic."
;
option
	author= "Carl Love <carl.j.love@gmail.com>, 20-May-2019",
	references= (                                                           #Ref #s:
		"https://en.wikipedia.org/wiki/Pearson%27s_chi-squared_test",         #*1
		"https://en.wikipedia.org/wiki/Yates%27s_correction_for_continuity",  #*2
		"https://en.wikipedia.org/wiki/G-test",                               #*3
		"https://en.wikipedia.org/wiki/Fisher%27s_exact_test",                #*4
		"Eric W Weisstein \"Fisher's Exact Test\" _MathWorld_"
		"--A Wolfram web resource:"
		"   http://mathworld.wolfram.com/FishersExactTest.html"               #*5
	)
;
uses AT= ArrayTools, St= Statistics, It= Iterator;
local
	C:= AT:-AddAlongDimension(O,1), R:= AT:-AddAlongDimension(O,2), #column & row sums
	c:= numelems(C), r:= numelems(R), #counts of columns & rows
	n:= add(R), #number of observations

	#All remaining code in `local` applies strictly to Fisher's.
	X, `&*!`:= mul@factorial~, 
	Cutoff:= &*!(O), #cut-off value for less-likely matrices
	CTs:= 0, #count of contingency tables with same row and column sums as O.
	Devs:= 0, #count of those tables with probabilities <= that of O.
	
	#Generate recursively all contingency tables whose row and column sums match O.
	#A very similar recursion is documented on the Maple help page
	# ?Iterator,BoundedComposition as example "Contingency Table".
	CountCTs:= (C,i)-> 
		`if`(i=r, 1, add(thisproc(C-X,i+1), X= It:-BoundedComposition(C,R[i]))),
		
        #Using the same recursion as above, generate the tables and compute their  
	#probabilities under the null hypothesis. The formula for the probability is
	#from reference *5. Sum probabilities of all those at least as improbable as O. 
	#Parameters: 
	#   C = column sums remaining to be filled; 
	#   F = product of factorials of entries of contingency table being built;
	#   i = row to be chosen this iteration
	AllCTs:= proc(C, F, i)
		if i = r then #Recursion ends.
			if irem(++CTs, 10^5) = 0 then Report() fi;
			#Last row is C, the unused portion of column sums:
			if (local P:= F*&*!(C)) >= Cutoff then ++Devs; 1/P else 0 fi
		else
			add(thisproc(C-X, F*&*!(X), i+1), X= It:-BoundedComposition(C, R[i]))
		fi
	end proc,
   
	Report:= ()-> userinfo(
		1, ContingencyTableTests, NoName, CTs, "tables evaluated, of which", Devs, 
		"are at least as improbable as the observed input."
	)     
;
	if method in {'Fisher', 'Fisher'['count']} then
		if method = 'Fisher'['count'] then return CountCTs() fi; 
		X:= AllCTs(C, 1, 1)*&*!(R)*&*!(C)/n!;  #see *5
		Report();
		X
	#Readers only interested in Fisher's exact test can ignore everything below
	#this point.
	else #Use one of the chi-squared methods: Pearson, G, Yates.
		local E::rtable; #matrix of expected values under null hypothesis 		
		try
			E:= rtable(
				(1..r, 1..c), 
				(i,j)-> R[i]*C[j],
				#The chi-squared tests aren't sufficiently
				#accurate if any cell has expected value < 5:
				datatype= satisfies(x-> x >= 5*n) 
			) / n;
		catch "unable to store":
			error
				"cell expected value too small for chosen method;"
				" use 'method'= 'Fisher'"
		end try;
		userinfo(
			1, ContingencyTableTests, 
			"Table of expected values under null hypothesis:", 
			print(
				`if`(
					E::'rtable'(posint), 
					E, 
					E=subsindets(E, Not(integer), evalf[1+ilog10(n)])
				)
      		)
		);
		#Pearson's, Yates's, and G all use a chi-square statistic,
		#each computed by a slightly different formula:
		local Chi2:= add@~table([
			'Pearson'= ((O,E)-> (O-E)^~2 /~ E),                     #see *1
			'Yates'= ((O,E)-> (abs~(O-E) -~ 1/2)^~2 /~ E),          #see *2
			'G'= ((O,E)-> 2*O*~map(x-> `if`(x=0, 0, ln(x)), O/~E))  #see *3
		]);
		1 - St:-CDF(ChiSquare((r-1)*(c-1)), Chi2[method](O,E)) 
	fi   
end proc:

 

[*3]Eric W. Weisstein, "Fisher's Exact Test", Mathworld--A Wolfram web resource.

You're trying to pack several matrices side-by-side into a single matrix, but the original matrices don't have the same number of rows, so it can't be done. Why not just export them individually?

Where you have k1 = 10, it should be k1:= 10.

It's done like this:

X:= Statistics:-RandomVariable(
   EmpiricalDistribution([-1, 0, 1, 2], probabilities= [0.2, 0.25, 0.35, 0.2])
):
Statistics:-Mean(X);

 

Your problem is simply a non-ASCII character in the last line of your code. Where there should be :=  you have :? where ? is some non-ASCII character that looks like :=

To plot it,

plots:-odeplot(dsn1, [w, x(w)])

and likewise for p(w).

The count for the probability's numerator is binomial(6,5)*binomial(9,5); for the denominator, binomial(15,5)*binomial(10,5)/3!.

I am posting this from my phone, where I have neither the ability to view your worksheet nor do my computation in Maple.

Assuming that your subscript notation just means another derivative, then it is easy to get these with implicitdiff. A single command can do all three:

EQ:= V^3-R*T*V^2/P-(B^2+P*B*R*T*sqrt(T)/(P-A))*V-A*B/(P*sqrt(T)) = 0:
indeps:= %op~([[T,V], [V,T], [V,V]]):
DP:=  table(value(indeps =~ map[3](%implicitdiff, EQ, P, indeps))); 

Now the three derivatives can be referenced by DP[T,V]DP[V,T], and DP[V,V]. You can very easily verify with simplify that the first two of these are equal.

First 143 144 145 146 147 148 149 Last Page 145 of 395