4
$\begingroup$

I'm trying to verify a solution to a simple probability problem using Mathematica. Here's the problem:

A drawer contains red socks and black socks. When two socks are drawn at random, the probability that both are red is 1/2. How small can the number of socks in the drawer be if the number of black socks is even?

From manual trial and error, it is trivial to find a solution of red = 15 and black = 6. However, when I try to verify the solution, Minimize and NMinimize fails to find an exact solution.

Minimize simply returns without having found a solution.

Minimize[{r + b, Binomial[r, 2]/Binomial[r + b, 2] == 1/2, r > 0, 
  b > 0, b/2 == i}, {r \[Element] Integers, b \[Element] Integers, 
  i \[Element] Integers}]

NMinimize throws an NMinimize:nosat error and proposes a close enough solution of red = 5 and black = 2. The probability works out to be $\frac{\binom{5}{2}}{\binom{7}{2}}=\frac{10}{21}$, which is not exactly 1/2.

NMinimize[{r + b, Binomial[r, 2]/Binomial[r + b, 2] == 1/2, r > 0, 
  b > 0, b/2 == i}, {r \[Element] Integers, b \[Element] Integers, 
  i \[Element] Integers}]

Is there something wrong with my setup, or is this a bug? I've tried all Methods, but neither made a difference

$\endgroup$

2 Answers 2

6
$\begingroup$

We solve the constrains at first.

sol = Solve[{Binomial[r, 2]/Binomial[r + b, 2] == 1/2, r > 0, b > 0, 
    b/2 == i}, {r ∈ Integers, b ∈ Integers, 
    i ∈ Integers}][[1]]

enter image description here

min=Minimize[{r + b /. sol /. C[1] -> c, 
   c ∈ PositiveIntegers}, {c}] // Simplify

{r, b} /. sol /. C[1] -> c /. min[[2]] // Simplify

{15, 6}

$\endgroup$
1
  • $\begingroup$ That works, thanks! What do you think makes Minimize struggle with this? Am I misusing it? $\endgroup$ Commented Sep 18, 2022 at 0:03
3
$\begingroup$

I suppose that the problem might be difficult for the minimizer because there are many local minima.

A direct search approach is the following:

MinimalBy[
  Tuples[Range[20], {2}],
  Apply[#1 + #2 + 1000 Boole@OddQ[#2] + 10000 Abs[(Binomial[#1, 2]/Binomial[#1 + #2, 2] - 1/2)]&]
]

(* Out: {{15, 6}} *)

One can get NMinimize to work without previously solving the constraints by modifying the objective function to minimize simultaneously the sum of r and b, and the absolute deviation from the desired probability. In practice, that means minimizing the product of the sum and the absolute value of (the probability - 1/2). Additionally, we provide all possible 2-tuples as starting points for the minimization method. Finally, the constraint that b should be even is treated directly through an auxiliary function rather than introduction another parameter.

Clear[even]
even[x_?NumericQ] := Boole@EvenQ[x]

NMinimize[{
  (r + b) Abs[Binomial[r, 2]/Binomial[r + b, 2] - 1/2],
  even[b] == 1,
  r \[Element] PositiveIntegers,
  b \[Element] PositiveIntegers},
  {r, b},
  Method -> {"NelderMead", "InitialPoints" -> Tuples[Range[20], {2}]}
]

(* Out: {0., {r -> 15, b -> 6}} *)
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.