- "markdown": "---\ntitle: \"Average Area of a Random Triangle\"\ndescription: \"An exercise in geometry, probability, and Python computing.\"\ndate: \"2025-03-10\"\ncategories: [probability, geometry, python]\nimage: \"image_2.png\"\n---\n\n\nThis post shows how to compute the expected value of the area of a random triangle analytically, numerically, and with Monte Carlo.\nSpecifically, suppose the vertices of a triangle are sampled uniformly from the perimeter of a circle with radius one (see picture).\nWhat is the expected value of the area of the triangle?\nIs it 0.5?\nIs it $\\pi$/5?\nIn this post, I'll share how I solved this problem.\n\n::: {#fig-triangles layout-ncol=2}\n\n\n\n\n\nTwo triangles inscribed in a circle. In our set-up, the blue circle has radius one and the vertices of the triangle are sampled uniformly on the perimeter of the circle. Our goal is to compute the average area of the triangle.\n:::\n\nFirst, I tried to find a closed-form expression for the expected value of the area.\nI expressed the $i$th triangle vertex as $(X_i,Y_i)$, where $X_i=\\cos(\\theta_i)$, $Y_i=\\sin(\\theta_i)$, and $\\theta_i\\sim \\textrm{Unif}(0,2\\pi)$ for $i=1,2,3$.\nI then searched on Google for the area of the triangle, and found the formula\n$$\nA=\\frac{1}{2}|X_1(Y_2-Y_3)+X_2(Y_3-Y_1)+X_3(Y_1-Y_2)|,\n$$\nwhich is half of the absolute value of the determinant of a matrix whose columns are $(X_1,X_2,X_3)$, $(Y_1,Y_2,Y_3)$, and $(1,1,1)$.\nNow, this formula depends on three angles, so I asked myself whether we could make a simplifying assumption that would reduce the dimension of the problem to two.\nAfter some thought, I realized that we can rotate the triangle so that $X_1=1$ and $Y_1=0$, which effectively sets $\\theta_1=0$ and leaves the area unchanged since rotation matrices have determinant one.\nThis simplifies the area formula to\n$$\nA=\\frac{1}{2}|Y_2-Y_3+X_2Y_3-X_3Y_2|.\n$$\nThe term $X_2Y_2-X_3Y_2$ looks like a trig identity, and so I searched on Google and found that\n$$\nX_2Y_3-X_3Y_2=-\\sin(\\theta_2-\\theta_3),\n$$\nallowing the area to be further simplified to\n$$\nA=\\frac{1}{2}|\\sin(\\theta_2)-\\sin(\\theta_3)-\\sin(\\theta_2-\\theta_3)|.\n$$\nI didn't see how to simplify this expression further, so I moved on to the expected value calculation.\n\nSince the angles $\\theta_2$ and $\\theta_3$ have joint density $1/(4\\pi^2)$, the expected value of the area is\n$$\nE[A]=\\frac{1}{8\\pi^2}\\int_0^{2\\pi}\\int_0^{2\\pi}|\\sin(\\theta_2)-\\sin(\\theta_3)-\\sin(\\theta_2-\\theta_3)|\\;d\\theta_2\\;d\\theta_3.\n$$\nAt this point I got stuck since I wasn't sure how to deal with the absolute value in the integrand.\nEventually I remembered that an absolute value can be written as a sum of a positive and a negative term, which applied to our problem leads to\n$$\nE[A]=\\frac{1}{8\\pi^2}\\left[\\int\\int_{f(\\theta_2,\\theta_3)>0}f(\\theta_2,\\theta_3)\\;d\\theta_2\\;d\\theta_3 - \\int\\int_{f(\\theta_2,\\theta_3)<0}f(\\theta_2,\\theta_3)\\;d\\theta_2\\;d\\theta_3 \\right],\n$$\nwhere\n$$\nf(\\theta_2,\\theta_3)=\\sin(\\theta_2)-\\sin(\\theta_3)-\\sin(\\theta_2-\\theta_3).\n$$\nThis looks a little better, and there's a hint of symmetry so I plotted the function and realized that 1) the second integral in the square brackets is the negative of the first integral and 2) the region over which $f(\\theta_2,\\theta_3)$ is positive is $\\theta_3>\\theta_2$.\nHere is a plot of the $f$ function, the code to generate the plot can be expanded by clicking on the arrow icon.\nNote that the axes limits are $(0,2\\pi)$.\n\n::: {#798acac6 .cell execution_count=1}\n``` {.python .cell-code code-fold=\"true\"}\nimport numpy as np\nimport matplotlib.pyplot as plt\n\n\ndef f(theta2: float, theta3: float) -> float:\n return np.sin(theta2) - np.sin(theta3) - np.sin(theta2 - theta3)\n\n\n# Create mesh grid\ntheta2 = np.linspace(0, 2 * np.pi, 100)\ntheta3 = np.linspace(0, 2 * np.pi, 100)\nTheta2, Theta3 = np.meshgrid(theta2, theta3)\nZ = f(Theta2, Theta3)\n\n# Plot contour\nplt.figure(figsize=(8, 6))\ncontour = plt.contourf(Theta2, Theta3, Z, levels=50, cmap=\"viridis\")\nplt.colorbar(contour)\nplt.xlabel(r\"$\\theta_2$\")\nplt.ylabel(r\"$\\theta_3$\")\nplt.title(\n r\"Contour Plot of $\\sin(\\theta_2) - \\sin(\\theta_3) - \\sin(\\theta_2 - \\theta_3)$\"\n)\nplt.show()\n```\n\n::: {.cell-output .cell-output-display}\n{}\n:::\n:::\n\n\nIt follows that the formula for the expected value of the area can be simplified to\n$$\nE[A]=\\frac{1}{4\\pi^2}\\int_0^{2\\pi}\\int_{\\theta_2}^{2\\pi}f(\\theta_2,\\theta_3)\\;d\\theta_3\\;d\\theta_2.\n$$\nI computed this last expression using WolframAlpha and found that the double integral equals $6\\pi$, which implies\n$$\n\\boxed{E[A]=\\frac{3}{2\\pi}}.\n$$\nWhat a pleasant surprise!\n\nWhile I was stuck on this integral, I computed the expected value numerically and with Monte Carlo.\nHere is the code:\n\n::: {#7d6d5198 .cell execution_count=2}\n``` {.python .cell-code}\nimport scipy.integrate as spi\n\n\ndef triangle_area(theta1: float, theta2: float) -> float:\n \"\"\"Compute the area of a triangle inscribed in a unit circle and with\n vertices at 0, theta1, and theta2 radians.\n \"\"\"\n X2 = np.cos(theta1)\n X3 = np.cos(theta2)\n Y2 = np.sin(theta1)\n Y3 = np.sin(theta2)\n area = 0.5 * np.abs(Y2 - Y3 + X2 * Y3 - X3 * Y2)\n return area\n\n\ndef compute_monte_carlo(n: int = 1_000_000, seed: int = 1) -> float:\n \"\"\"Compute expected value of triangle area using Monte Carlo.\"\"\"\n np.random.seed(seed)\n theta1, theta2 = np.random.uniform(low=0, high=2 * np.pi, size=(2, n))\n areas = triangle_area(theta1, theta2)\n return np.mean(areas)\n\n\ndef compute_numerical_integral() -> float:\n \"\"\"Compute expected value of triangle area using numerical integration.\"\"\"\n result, _ = spi.dblquad(\n triangle_area, 0, 2 * np.pi, lambda x: 0, lambda x: 2 * np.pi\n )\n result /= 4 * np.pi**2\n return result\n\n\n# Compute expected value using both methods\nintegral_est = compute_numerical_integral()\nmc_est = compute_monte_carlo()\n```\n:::\n\n\nLet's put the results in a table to compare them.\n\n::: {#tbl-expected-area .cell tbl-align='center' tbl-cap='Comparison of values for the average area.' execution_count=3}\n``` {.python .cell-code code-fold=\"true\"}\nimport pandas as pd\nfrom IPython.display import Markdown\n\ndata = {\n \"Method\": [\"Exact Value\", \"Numerical Integration\", \"Monte Carlo\"],\n \"Expected Area\": [3 / 2 / np.pi, integral_est, mc_est],\n}\ndf = pd.DataFrame(data)\n\n# Format numeric values to three significant figures\ndf = df.map(lambda x: f\"{x:.4g}\" if isinstance(x, float) else x)\n\nmarkdown_table = df.to_markdown(index=False)\n\nstyled_table = f\"\"\"\n<div style=\"display: inline-block; width: 40%; min-width: 200px;\">\n{markdown_table}\n </div>\n\"\"\"\n\ndisplay(Markdown(styled_table))\n```\n\n::: {.cell-output .cell-output-display .cell-output-markdown}\n\n<div style=\"display: inline-block; width: 40%; min-width: 200px;\">\n| Method | Expected Area |\n|:----------------------|----------------:|\n| Exact Value | 0.4775 |\n| Numerical Integration | 0.4775 |\n| Monte Carlo | 0.4773 |\n </div>\n\n:::\n:::\n\n\nIt's great to see all three methods give the same answer.\nOverall, this was a fun exercise and gave me an excuse to get some experience creating a blog with Quarto!\n\n",
0 commit comments