2
$\begingroup$

I draw the attached plot drawn with the following code. Each color represents a cell. How can I find the average number of cells that one user (any point within the square) can be in? By average I mean "the user can be anywhere within the whole area". For example, averaging over 10000 random locations of a point within the square area.

enter image description here

 ell = t [Function] t^-20.75];500, 2}]; 
 H = RandomVariate[ExponentialDistribution[1], 50]; 
 r = Table[H[[i]] ell[Norm[{x, y} - X[[i]]]], {i, 1, 20}]; 
 s = Table[r[[i]]/( Total[Delete[r, i]] + 30.99), {i, 1, 20}]; 

 Show[Table[
 RegionPlot[s[[i]] >= 0.1, {x, -0.7, 0.7}, {y, -0.7, 0.7}, 
 MaxRecursion -> 2, PlotPoints -> 40, 
 PlotStyle -> {FaceForm[{Opacity[0.5], ColorData[97][i]}]}, 
 BoundaryStyle -> ColorData[97][i]], {i, 1, 20}], Graphics[Point[X]]]
Improve this question
$\endgroup$
10
  • $\begingroup$ Wouldn't it be easier (and more correct) to sum up the areas of the shaded regions (with multiplicity) and to divide the result by the area of the square? $\endgroup$
    – Henrik Schumacher
    Jul 30, 2020 at 10:27
  • $\begingroup$ @HenrikSchumacher, seems to be a very good idea. But I do not know how to do it. I really appreciate if you can help me on this. $\endgroup$
    – MGK
    Jul 30, 2020 at 10:39
  • 1
    $\begingroup$ Sum[Area@ImplicitRegion[f >= 0.1, {{x, -0.7, 0.7}, {y, -0.7, 0.7}}],{f, s}]/(2 0.7)^2? See also the documentation of Area for a detailed explanation of the options of Area that allow you to control the accuracy of the computation (AccuracyGoal and PrecisionGoal). $\endgroup$
    – Henrik Schumacher
    Jul 30, 2020 at 10:45
  • $\begingroup$ @HenrikSchumacher, Thanks a lot. It seems to be computationally very expensive approach. Is there any alternative way that you can think of!! $\endgroup$
    – MGK
    Jul 30, 2020 at 11:14
  • 1
    $\begingroup$ Why do you think so? It is not more expensive to compute an approximate area of an impicit region than computing an approximate triangle mesh via RegionPlot. Basically that is what Area does: It meshes the domain and sums the triangle areas. $\endgroup$
    – Henrik Schumacher
    Jul 30, 2020 at 11:17

1 Answer 1

Reset to default
2
$\begingroup$

Well, a Monte-Carlo approach could look like this. There is a considerable amount of number crunching involved, so we better compile the working horse function:

cf = Compile[{{H, _Real, 1}, {X1, _Real, 1}, {X2, _Real, 1}, {Z, _Real, 1}},
   Block[{n, R, r, x1, x2},
    x1 = Compile`GetElement[Z, 1];
    x2 = Compile`GetElement[Z, 2];
    r = H ((X1 - x1)^2 + (X2 - x2)^2)^-1.875;
    R = Total[r] + 3.99;
    Total[UnitStep[r - (0.1/1.1) R]]
    ]
   ,
   CompilationTarget -> "C",
   RuntimeAttributes -> {Listable},
   Parallelization -> True,
   RuntimeOptions -> "Speed"
   ];

Now:

pts = RandomReal[{-0.7, 0.7}, {10000000, 2}];
{X1, X2} = Transpose[X];
N@Mean@cf[H, X1, X2, pts]

2.05869

The function cf counts how many entries of the vector s are above 0.1 for a given {x,y} = Z. Applying cf in a listable way to all elements of pts and taking the mean yields the empirical average number of satisfied inequalities.

Improve this answer
$\endgroup$

Your Answer

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

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