Carl Love

Carl Love

28055 Reputation

25 Badges

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

MaplePrimes Activity


These are replies submitted by Carl Love

@Christian Wolinski So, you managed to put the OP's code into an executable format, but you didn't think to upload it so that others could continue the work without needing to start from scratch? Collaborate.

@asma khan If matrices A and B have the same number of rows, then they can be attached side-by-side as Tom showed:

<A|B>

If they have the same number of columns, then they can be attached top-to-bottom by

<A,B>

These commands are not limited to two matrices. They can be interchangeably nested to arbitrary depth to create matrices of arbitrary complexity.

@nm I do not "object" to you saying "internal" error. I merely thought that you'd misinterpretted the int in `latex/int` to stand for "internal". On the other hand, with your clarification on your definition of "internal error", it seems to mean the same thing as "bug".

As usual, I'm just trying to offer a practical solution so that a user can proceed with their work. I'm not trying to excuse the error. I'm not addressing what should be; I'm addressing what is. I do think, like you, that the command "should  internally handle these issues", but I don't work for Maplesoft! I don't have the power to fix bugs! You should file an SCR. You and I have had this conversation several times before. Please do not respond to my workarounds as if I had the power to make things work as they "should". It is very annoying to me when you do that. I am just a "user" who likes answering questions and providing help.

Practical means, among other things, that the solution may not work on 10,000 equations all at once without error. I've never heard of anyone performing a large project with the latex command without needing to do some "hand tuning", and I'm amazed that you've not had a problem with the command before. It's a very weak feature of Maple, IMO.

If you'd like to correct the bug before someone at Maplesoft gets to it, the procedure `latex/int` would be trivial to correct by taking out the first two lines, which check that the integral has exactly 2 arguments.

@David Sycamore Got it. It wasn't clear from your original whether primality or parity should be checked first (and it only makes a difference for n=2). Also, you wrote a(n+1) = n+2 where it should be a(n+1) = a(n)+2. So with that in mind, the procedure can be simplified as follows. I also show one way you can get the pairwise output that you want. There are many ways that this could be formatted.

restart:
A:= proc(n::posint) 
option remember; 
local a:= thisproc(n-1)+1;   
   a +`if`(a::prime, 0, 2-irem(a,2))
end proc:
A(1):= 1: 
:
seq('A'(n)=A(n), n= 1..20);
   A(1) = 1, A(2) = 2, A(3) = 3, A(4) = 6, A(5) = 7, A(6) = 10, 
   A(7) = 11, A(8) = 14, A(9) = 16, A(10) = 17, A(11) = 20, 
   A(12) = 22, A(13) = 23, A(14) = 26, A(15) = 28, A(16) = 29, 
   A(17) = 32, A(18) = 34, A(19) = 36, A(20) = 37

And there's no reason that such a list must start at n=1. You could just as well do n= 100..120.

@mmcdara Yes, I'd say that it's never wise to use an indexed name such as x[0] together with its parent x when a subscripted name such as x__0 will do the job. The main problem arises when assignments are made. Any assignment made to x[...when x is unassigned will turn into a table. Since tables have the last name evaluation property, it'll still be possible to some extent to use x as a variable. Procedures also use last name evaluation, which is why the OP's use of Zeta doesn't lead to total nonsense. On the other hand, an unwary assignment to x will cause x[...] to be nonsense. 

I'm not on my computer right now, just on my phone, and I don't know much about such problems anyway. Yet your r < a sticks out to me. Shouldn't that be r <= 1, if anything at all? And how about r >= 0?

@mmcdara Your solution works for expanded polynomials (which is the case for this Question), but fails for almost anything more complicated such as rational functions or transcendental functions. Besides, there is a specific command to handle this precise situation gracefully: evalc.

@Rouben Rostamian  Upon carefully rereading the Question, I think that your suggestion fulfills the OP's intention. I at first thought that they wanted to integrate the right side and maintain it as an equation with the same left side, which is what my code does.

@vv Nothing has been deleted. If you click on it from the main menus, you may see Page Not Found. This is just a bug in MaplePrimes. If instead you go to the "Active Conversations" page, then click on it, you will see it. It seems to take MaplePrimes somewhere between 1 and 2 hours to re-index the main menu links. After that, it'll appear as normally.

@vv Although it wasn't I who branched this off from its former parent Question, I think that it was the right decision to do so:

  • It's an explicit mathematical question.
  • Maple can help with its solution.
  • An elaborate Answer using Maple has been posted.

@Rouben Rostamian To do that, you'd need to implement the arithmetic operations multiplication, addition, and greater than in Logic. It would be quite elaborate, but in theory doable[*1]. If a problem can be solved with ordinary arithmetic, it'd probably be best to avoid using Logic.

Note that the only arithmetic operations needed for Joe's solution are less than, add 1, and subtract 1.

[*1]as shown with great detail, precision, and fanfare in Whitehead and Russell's Principia Mathematica, 1913.

@q41b3 You're welcome. I don't know if it's worth studying your book for its Maple content due to the age. Perhaps only the math content is relevant anymore. There's been a lot added to Maple since 2002, including the packages LinearAlgebraStatistics, and Finance (which handles the main stochastic differential equations used in finance).

I've extensively updated this code with comments and to use syntax features new to Maple 2019. In particular, I like the new local syntax, which allows placement of the declaration close to where the variable is used, which I think is more robust and easier to read. I've also used the new increment syntax, ++var, and the new neutral-operator syntax, which allows the use of the operator in functional form (i.e., coming before it arguments) without the use of back quotes.

I recently encountered on MaplePrimes a problem for which Fisher's exact test is ideal, and I've included that problem.

Since MaplePrimes doesn't reproduce Code Edit Regions, here's my code, and the worksheet follows.

#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'}
	#Use 'method' = 'Fisher[count]' to count how many table variations are 
	#needed for Fisher's.
)
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>, 29-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.
	`&*!`:= 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:= proc(C, i)
	local X; 
		`if`(i=r, 1, add(thisproc(C-X,i+1), X= It:-BoundedComposition(C,R[i])))
	end proc,
		
        #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:= (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
			local X;
			add(thisproc(C-X, F*&*!(X), i+1), X= It:-BoundedComposition(C,R[i]))
		fi,
   
	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(C,1) fi; 
		local p:= AllCTs(C, 1, 1)*&*!(R)*&*!(C)/n!;  #p-value (see *5)
		Report();
		p
	#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]/n,
				#The chi-squared tests aren't sufficiently
				#accurate if any cell has expected value < 5:
				datatype= satisfies(x-> x >= 5) 
			)
		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)) #p-value 
	fi   
end proc:


 

Contingency table tests, including Fisher's exact test

Author: Carl Love <carl.j.love@gmail.com>, 29-May-2019

restart:

This code uses features new to Maple 2019.

Problem: On an episode of the television show Survivor, there were 15 contestants divided into 2 teams of sizes 9 and 6. The producers decided to divide them (using an unknown selection process) into 3 teams of size 5. The resulting selection was that one of the new teams was drawn entirely from the size-9 team, a second was drawn entirely from the size-6 team, and, of course, the third was comprised of the remaining 5 contestants. Assuming that the selection process was random (this being the null hypothesis), what is the probability of achieving a result at least this extreme?

 

Solution: Although this probability can be easily computed by hand using basic combinatorics (see the end of this worksheet), it's also a good problem for Fisher's exact test. In the matrix below, the rows show the division into the original teams, and the columns show the division into the new teams.

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

infolevel[ContingencyTableTests]:= 1:

ContingencyTableTests(SurvivorTeams, method= Fisher): % = evalf(%);

25 tables evaluated, of which 6 are at least as improbable as the observed input.

6/1001 = 0.5994005994e-2

Since the number of tables which need to be evaluated can reach into the millions even for modestly sized problems, I've included an option to count them prior to doing the full computation:

ContingencyTableTests(SurvivorTeams, method= Fisher[count]);

25

This problem is too small to use any of the more-standard contingency table tests based on a chi-squared statistic, and my program can detect that:

ContingencyTableTests(SurvivorTeams);

Error, (in ContingencyTableTests) cell expected value too small for chosen method; use 'method'= 'Fisher'

There are numerous intuitive ways to do the computation by hand. Here's one:

binomial(9,5)*binomial(6,5)/(binomial(15,5)*binomial(10,5)/3!);

6/1001

 


 

Download FishersExact.mw

@q41b3 Well, you can still choose this Answer as best.

@q41b3 The difference between else if and elif is purely syntactic. Each if must have a corresponding terminator (fiend if, or end---those three being equivalent in this context). Using else if starts a new if statement, so it requires two terminators. 

First 271 272 273 274 275 276 277 Last Page 273 of 709