Monte Carlo Decision Analysis Engine
Build a pure Python Monte Carlo engine to optimize investment allocations.
Introduction
30 Second Summary
Every major investment decision comes down to one question: what is the best way to allocate resources when the future is uncertain? This is the exact problem that quantitative analysts solve with computational simulation.
In this project, you will build a Monte Carlo Decision Analysis Engine entirely from Python's standard library. Your engine will model probability distributions, run thousands of simulated investment trials, and identify the optimal way to allocate a $10M fund across three strategies with uncertain returns.
What You'll Build
A single-file Python analysis engine that solves a complex portfolio allocation problem under uncertainty, exercising 30+ distinct language constructs as thinking tools for analytical problem-solving.
By the end of this project, you'll have:
- A distribution modeling layer using abstract base classes, frozen dataclasses, operator overloading, and singledispatch factories.
- A simulation engine powered by generators, itertools, custom decorators, and context managers.
- A decision analysis framework using structural pattern matching and combinatoric search to find optimal allocations.
- A reporting system with statistical analysis, correlation matrices, and ASCII histogram visualization.
- Secret Mission: Add sensitivity analysis using linear regression and implement a risk parity optimizer that balances risk contributions across strategies.
Are there any prerequisites?
You should be comfortable with Python fundamentals: functions, classes, loops, and basic data structures. Familiarity with basic probability and statistics concepts (mean, standard deviation) will also help you follow along.
Before We Start
You are about to build a Monte Carlo Decision Analysis Engine that helps a venture capital firm allocate a $10M fund across three investment strategies with uncertain returns.
The engine models probability distributions, runs thousands of simulated trials, and finds the optimal allocation that maximizes risk-adjusted returns. You will build this using nothing but Python's standard library, which forces you to exercise the language's full power.
In this step, get ready to:
- Verify you have Python 3.12+ installed.
- Create your project directory and open it in VS Code.
- Create the monte_carlo_engine.py file.
Verify your Python version
This project uses Python language features introduced in version 3.12, including the type statement syntax for type aliases. You need Python 3.12 or higher to run the engine.
- Open your terminal by pressing Cmd+Space (macOS) or the Windows key (Windows), typing Terminal (macOS) or cmd (Windows), and pressing Enter.
- Check your Python version by running this command:
python3 --version✔️ I see version 3.12 or higher
You should see output like Python 3.12.x or Python 3.13.x or Python 3.14.x.
Your Python version is ready. Move on to the next substep.
ⓧ I see an older version
If you see Python 3.11 or lower, you need to install a newer version. On macOS, the system python3 may point to an older version bundled with Apple's developer tools.
- Go to python.org/downloads and download the latest Python installer for your platform.
- On macOS, download the macOS 64-bit universal2 installer and run it.
- On Windows, download the Python install manager and run it.
- After installing, close your terminal and open a new one.
- Run python3 --version again to confirm you see 3.12 or higher.
Why does macOS show an older version?
macOS ships with an Apple-managed python3 in /usr/bin/python3 that is typically older. Installing from python.org places a newer version in /usr/local/bin/ which takes priority on your PATH after a terminal restart.
ⓧ Command not found
If you see an error about the command not being found, Python is not installed on your system yet.
- Go to python.org/downloads and download the latest Python installer for your platform.
- On macOS, download the macOS 64-bit universal2 installer, double-click to launch it, and follow the prompts.
- On Windows, download the Python install manager and run it.
- After installing, close your terminal and open a new one.
- Run python3 --version again to confirm you see 3.12 or higher.
Still seeing command not found after installing?
On Windows, try running python --version (without the 3) instead. The Windows installer may register the command as python rather than python3.
If neither works, your terminal may need a full restart. Close all terminal windows and reopen.
Help me troubleshoot my Python installation.
Create your project directory and file
- Move to your Desktop by running this command:
cd ~/Desktop- Create your project folder and move into it by running these commands:
mkdir monte-carlo-engine
cd monte-carlo-engine- Create the main Python file by running:
touch monte_carlo_engine.py- Confirm the file exists by running:
lsYou should see monte_carlo_engine.py in the output.
- Open the project folder in VS Code by pressing Cmd+Space (macOS) or the Windows key (Windows), typing VS Code, and pressing Enter.
- In VS Code, click File in the top menu bar.
- Click Open Folder.
- Navigate to your Desktop, select the monte-carlo-engine folder, and click Open.
- In the Explorer sidebar on the left, click monte_carlo_engine.py to open it in the editor.
Why a single file?
This project lives entirely in one file. The goal is to focus on Python's language features rather than project packaging or multi-module imports. You will build up the engine section by section, and each section demonstrates a different set of Python constructs.
Don't see the file in VS Code?
Make sure you opened the monte-carlo-engine folder (not your entire Desktop). The file monte_carlo_engine.py should appear in the Explorer sidebar.
If the sidebar is hidden, press Cmd+B (macOS) or Ctrl+B (Windows) to toggle it.
Help me find my project file in VS Code.
Your environment is ready. Next up, you will build the distribution modeling layer that represents uncertain investment returns using abstract classes, frozen dataclasses, and operator overloading.
Set Up Imports and Project File
You're building a Monte Carlo Decision Analysis Engine using nothing but Python's standard library. Every module you import earns its place by powering a specific layer of the engine.
In this step, you'll add all the import statements that bring in probability sampling, statistical analysis, combinatoric search, and architectural patterns. Python's standard library covers all of it.
In this step, get ready to:
- Add the full set of standard library imports that power the simulation engine.
- Enable forward-reference annotations for cleaner type hints.
Add the standard library imports
The engine pulls from 11 different standard library modules. Each one handles a specific concern, from random sampling to functional programming patterns to abstract class definitions.
- In monte_carlo_engine.py, add the following imports:
from __future__ import annotations
import random
import statistics
import itertools
import functools
import contextlib
import collections
import enum
import math
import time
from dataclasses import dataclass, field
from typing import Protocol, TypeVar, Iterator, Callable
from abc import ABC, abstractmethodWhat does this code do?
from __future__ import annotations enables postponed evaluation of type annotations. This lets you reference classes before they're defined (forward references).
The imports divide into four groups based on what they power in the engine:
- Simulation: random for probability sampling, statistics for analysis, math for numerical operations, time for performance timing.
- Architecture: itertools for combinatoric iteration, functools for decorators and caching, contextlib for context managers.
- Data modeling: dataclasses for structured data, collections for specialized containers, enum for named constants.
- Type system: typing for Protocol and generics, abc for abstract base classes.
- Save monte_carlo_engine.py.
- Confirm there are no red underlines or error squiggles in VS Code. All these modules are built into Python's standard library, so no installation is needed.
Seeing red underlines on an import?
- Check that VS Code is using your Python 3.12+ interpreter. Click the Python version in the bottom status bar and select the correct one.
- Make sure from __future__ import annotations is the very first line in the file. Nothing can come before a __future__ import except comments.
Help me fix import errors in my Python file.
✔️ Awesome, I've got everything!
Great. Double check you've saved your file.
ⓧ I'd like to double check the full code
Here's exactly what monte_carlo_engine.py should contain at this point:
from __future__ import annotations
import random
import statistics
import itertools
import functools
import contextlib
import collections
import enum
import math
import time
from dataclasses import dataclass, field
from typing import Protocol, TypeVar, Iterator, Callable
from abc import ABC, abstractmethodWhy standard library only?
By avoiding external libraries like NumPy or pandas, every Python language feature (generators, decorators, pattern matching, operator overloading) earns its place by solving a real problem. There's no shortcut API to lean on.
This forces you to use Python's full syntax as a thinking tool, not just as glue between library calls.
Your project file is loaded with everything you need. Next up, you'll build the distribution modeling layer that gives the engine its probability-powered foundation.
Build the Distribution Modeling Layer
Your imports are in place. Now you need a way to model uncertain investment returns.
In this step, you'll build the foundation of the engine: a distribution modeling layer that uses abstract classes, frozen dataclasses, operator overloading, and a factory pattern. This layer gives the simulation engine a clean interface for generating random samples from different probability distributions.
In this step, get ready to:
- Define custom exceptions, an enum, a protocol, and type aliases for the type system.
- Create an abstract base class and three frozen dataclass distributions.
- Add operator overloading and a singledispatch factory for flexible distribution creation.
Define exceptions, enum, protocol, and type aliases
Every modeling layer needs a clear error hierarchy and type definitions. Custom exceptions give precise feedback when a distribution receives invalid parameters. An enum encodes risk levels with associated multipliers. A protocol defines the structural contract that all distributions must satisfy.
- At the very top of monte_carlo_engine.py, above from __future__ import annotations, add this module docstring:
"""\nMonte Carlo Decision Analysis Engine\nPure Python Standard Library Implementation\n"""- Below all your import statements, add the Part 1 section header and custom exception classes by copying this code:
\n# ============================================================\n# PART 1: Distribution Modeling Layer\n# ============================================================\n\n# --- Custom Exceptions ---\n\nclass DistributionError(Exception):\n """Base exception for distribution errors."""\n pass\n\n\nclass InvalidParameterError(DistributionError):\n """Raised when distribution parameters are invalid."""\n\n def __init__(self, param: str, value: float, reason: str) -> None:\n self.param = param\n self.value = value\n super().__init__(f"Invalid {param}={value}: {reason}")What does this code do?
DistributionError is the base exception for anything that goes wrong with distributions. InvalidParameterError inherits from it and stores the parameter name, value, and reason so error messages are precise.
The __init__ override calls super().__init__() with a formatted string, giving you messages like "Invalid sigma=−0.5: must be positive" when something breaks.
- Save monte_carlo_engine.py.
- You won't see output yet. The next chunk adds the risk level enum.
- Below the exception classes, add the RiskLevel enum by copying this code:
\n# --- Enums ---\n\nclass RiskLevel(enum.Enum):\n LOW = "low"\n MEDIUM = "medium"\n HIGH = "high"\n EXTREME = "extreme"\n\n @property\n def multiplier(self) -> float:\n multipliers = {\n RiskLevel.LOW: 0.5,\n RiskLevel.MEDIUM: 1.0,\n RiskLevel.HIGH: 1.5,\n RiskLevel.EXTREME: 2.0,\n }\n return multipliers[self]What does this code do?
RiskLevel is a Python Enum with four members. Each member maps to a string value like "low" or "extreme".
The @property decorator turns multiplier into a computed attribute. When you access RiskLevel.EXTREME.multiplier, it returns 2.0 without needing parentheses.
- Save monte_carlo_engine.py.
- You won't see output yet. The next chunk defines the protocol and type aliases.
- Below the RiskLevel enum, add the protocol definition, type variable, and type aliases by copying this code:
\n# --- Protocol (structural typing) ---\n\nT = TypeVar("T", bound="Distribution")\n\n\nclass Samplable(Protocol):\n def sample(self, n: int, *, seed: int | None = None) -> Iterator[float]: ...\n\n\n# --- Type Alias (3.12+ syntax) ---\n\ntype Samples = list[float]\ntype Parameters = dict[str, float]What does this code do?
TypeVar creates a generic type variable T bounded to Distribution. This lets type checkers understand that methods returning T will return the same subclass they were called on.
The Samplable Protocol defines structural typing. Any class with a matching sample() signature satisfies the protocol without explicitly inheriting from it.
The type statement is Python 3.12+ syntax for creating type aliases. It replaces the older TypeAlias annotation pattern with a cleaner declaration.
- Save monte_carlo_engine.py.
- You won't see output yet. The next substep builds the abstract base class and concrete distributions.
Seeing a syntax error on the type statement?
The type Samples = list[float] syntax requires Python 3.12 or higher. Run python3 --version to confirm your version. If you're on 3.11 or lower, you need to upgrade.
help me fix a type alias syntax error in Python
Create the Distribution ABC and frozen dataclass distributions
The abstract base class defines the interface that every distribution must implement. Frozen dataclasses then provide immutable, type-safe implementations. The frozen=True flag prevents mutation after creation, which is important for reproducible simulations.
- Below the type aliases, add the Distribution abstract base class by copying this code:
\n# --- Abstract Base Class ---\n\nclass Distribution(ABC):\n """Abstract base for probability distributions."""\n\n @abstractmethod\n def sample_one(self, rng: random.Random) -> float:\n """Generate a single sample."""\n ...\n\n @abstractmethod\n def expected_value(self) -> float:\n """Theoretical mean."""\n ...\n\n @abstractmethod\n def describe(self) -> str:\n """Human-readable description."""\n ...\n\n def sample(self, n: int, *, seed: int | None = None) -> Iterator[float]:\n """Generate n samples as a lazy iterator."""\n rng = random.Random(seed)\n for _ in itertools.repeat(None, n):\n yield self.sample_one(rng)\n\n def __repr__(self) -> str:\n return self.describe()What does this code do?
The ABC (Abstract Base Class) forces every subclass to implement sample_one(), expected_value(), and describe(). If a subclass forgets one, Python raises TypeError at instantiation.
The concrete sample() method is a generator. It uses itertools.repeat(None, n) instead of range(n) to loop n times. The yield keyword makes this a lazy iterator that produces samples one at a time without storing them all in memory.
- Save monte_carlo_engine.py.
- You won't see output yet. The next chunk implements the first concrete distribution.
- Below the Distribution class, add the frozen dataclass section header and the NormalReturn distribution by copying this code:
\n# --- Frozen Dataclass Distributions ---\n\n@dataclass(frozen=True, slots=True)\nclass NormalReturn(Distribution):\n """Normally distributed returns."""\n\n mu: float\n sigma: float\n label: str = field(default="Normal", repr=False)\n\n def __post_init__(self) -> None:\n if self.sigma <= 0:\n raise InvalidParameterError("sigma", self.sigma, "must be positive")\n\n def sample_one(self, rng: random.Random) -> float:\n return rng.gauss(self.mu, self.sigma)\n\n def expected_value(self) -> float:\n return self.mu\n\n def describe(self) -> str:\n return f"{self.label}(\u03bc={self.mu:.2%}, \u03c3={self.sigma:.2%})"What does this code do?
The @dataclass(frozen=True, slots=True) decorator generates __init__, __eq__, and __hash__ automatically. frozen=True makes instances immutable. slots=True uses __slots__ instead of a __dict__ for faster attribute access and lower memory usage.
__post_init__ runs after the generated __init__ completes. It validates that sigma is positive, raising the custom exception immediately if someone creates NormalReturn(0.1, -0.5).
The field(default="Normal", repr=False) on label gives it a default value and hides it from the auto-generated __repr__.
- Save monte_carlo_engine.py.
- You won't see output yet. The next chunk adds the triangular and log-normal distributions.
- Below the NormalReturn class, add the TriangularReturn and LogNormalReturn distributions by copying this code:
\n@dataclass(frozen=True, slots=True)\nclass TriangularReturn(Distribution):\n """Triangular distribution for bounded uncertain returns."""\n\n low: float\n mode: float\n high: float\n label: str = field(default="Triangular", repr=False)\n\n def __post_init__(self) -> None:\n if not (self.low <= self.mode <= self.high):\n raise InvalidParameterError(\n "mode", self.mode,\n f"must satisfy low({self.low}) <= mode <= high({self.high})"\n )\n\n def sample_one(self, rng: random.Random) -> float:\n return rng.triangular(self.low, self.high, self.mode)\n\n def expected_value(self) -> float:\n return (self.low + self.mode + self.high) / 3\n\n def describe(self) -> str:\n return f"{self.label}(low={self.low:.2%}, mode={self.mode:.2%}, high={self.high:.2%})"\n\n\n@dataclass(frozen=True, slots=True)\nclass LogNormalReturn(Distribution):\n """Log-normal distribution for multiplicative returns."""\n\n mu: float\n sigma: float\n label: str = field(default="LogNormal", repr=False)\n\n def __post_init__(self) -> None:\n if self.sigma <= 0:\n raise InvalidParameterError("sigma", self.sigma, "must be positive")\n\n def sample_one(self, rng: random.Random) -> float:\n return math.exp(rng.gauss(self.mu, self.sigma)) - 1\n\n def expected_value(self) -> float:\n return math.exp(self.mu + self.sigma**2 / 2) - 1\n\n def describe(self) -> str:\n return f"{self.label}(\u03bc={self.mu:.3f}, \u03c3={self.sigma:.3f})"What does this code do?
TriangularReturn models returns bounded between a low and high with a most-likely mode. Its __post_init__ validates that low <= mode <= high. The expected value is the simple average of all three parameters.
LogNormalReturn models multiplicative growth (like compound returns). It samples from a normal distribution in log-space using math.exp(rng.gauss(...)) - 1, which naturally produces the right-skewed returns you see in venture capital outcomes.
- Save monte_carlo_engine.py.
- You won't see output yet. The next substep adds operator overloading and a factory function.
Seeing a TypeError about frozen classes?
If Python raises FrozenInstanceError inside __post_init__, check that you are only reading self.sigma and self.mode for validation (not assigning to them). Frozen dataclasses prevent attribute assignment after init, but __post_init__ can still read attributes.
help me fix a FrozenInstanceError in my dataclass __post_init__
Add operator overloading and the singledispatch factory
Operator overloading lets you combine distributions using + and * syntax. A singledispatch factory provides a single entry point that creates the right distribution type based on the input format.
- In monte_carlo_engine.py, scroll to the NormalReturn class and find the describe method.
- Find this line:
def describe(self) -> str:
return f"{self.label}(\u03bc={self.mu:.2%}, \u03c3={self.sigma:.2%})"- Directly below this describe method (still inside the NormalReturn class), add these three operator overloading methods:
def __add__(self, other: NormalReturn) -> NormalReturn:\n """Sum of independent normals is normal."""\n combined_mu = self.mu + other.mu\n combined_sigma = math.sqrt(self.sigma**2 + other.sigma**2)\n return NormalReturn(combined_mu, combined_sigma, f"{self.label}+{other.label}")\n\n def __mul__(self, scalar: float) -> NormalReturn:\n """Scale a distribution."""\n return NormalReturn(self.mu * scalar, self.sigma * abs(scalar), self.label)\n\n def __rmul__(self, scalar: float) -> NormalReturn:\n return self.__mul__(scalar)What does this code do?
__add__ lets you write dist_a + dist_b to combine two independent normal distributions. The sum of independent normals is normal, with combined mean and combined variance (via Pythagorean theorem on standard deviations).
__mul__ lets you write dist * 2.0 to scale a distribution. __rmul__ handles the reversed case 2.0 * dist. Without __rmul__, Python wouldn't know how to multiply a float by your custom class.
- Save monte_carlo_engine.py.
- You won't see output yet. The next chunk adds the singledispatch factory.
- Below the LogNormalReturn class (after all three distribution classes), add the singledispatch factory by copying this code:
\n# --- @singledispatch factory ---\n\n@functools.singledispatch\ndef create_distribution(spec) -> Distribution:\n """Factory: create distribution from various input types."""\n raise TypeError(f"Cannot create distribution from {type(spec).__name__}")\n\n\n@create_distribution.register(dict)\ndef _from_dict(spec: dict) -> Distribution:\n match spec:\n case {"type": "normal", "mu": mu, "sigma": sigma, **rest}:\n return NormalReturn(mu, sigma, rest.get("label", "Normal"))\n case {"type": "triangular", "low": low, "mode": mode, "high": high, **rest}:\n return TriangularReturn(low, mode, high, rest.get("label", "Triangular"))\n case {"type": "lognormal", "mu": mu, "sigma": sigma, **rest}:\n return LogNormalReturn(mu, sigma, rest.get("label", "LogNormal"))\n case _:\n raise DistributionError(f"Unknown distribution spec: {spec}")\n\n\n@create_distribution.register(tuple)\ndef _from_tuple(spec: tuple) -> Distribution:\n """Quick creation: (low, mode, high) -> Triangular."""\n if len(spec) == 3:\n return TriangularReturn(*spec)\n raise DistributionError(f"Tuple must have 3 elements, got {len(spec)}")What does this code do?
@functools.singledispatch creates a generic function that dispatches to different implementations based on the type of the first argument. Call create_distribution({"type": "normal", ...}) and Python routes to _from_dict. Call create_distribution((-0.05, 0.08, 0.25)) and it routes to _from_tuple.
Inside _from_dict, the match/case statement uses mapping patterns to destructure the dictionary. The **rest capture grabs any extra keys (like label). The wildcard case _ handles unrecognized specs.
Why singledispatch instead of if/elif?
Singledispatch is extensible. If you add a new input format later, you register a new handler without modifying the original function. It follows the open/closed principle: open for extension, closed for modification.
- Save monte_carlo_engine.py.
- Run a syntax check by running this command in your terminal:
python3 monte_carlo_engine.pyYou should see no output and no errors. The file defines classes and functions but has no main() call yet, so Python parses everything and exits silently. This confirms there are no syntax errors.
Seeing a SyntaxError on match/case?
Structural pattern matching (match/case) requires Python 3.10+ and the type alias syntax requires 3.12+. Run python3 --version and confirm you see 3.12 or higher.
help me resolve a SyntaxError with match/case in Python
✔️ Awesome, I've got everything!
Great. Double check you've saved monte_carlo_engine.py and the syntax check passed with no output.
ⓧ I'd like to double check the full code
"""\nMonte Carlo Decision Analysis Engine\nPure Python Standard Library Implementation\n"""\n\nfrom __future__ import annotations\n\nimport random\nimport statistics\nimport itertools\nimport functools\nimport contextlib\nimport collections\nimport enum\nimport math\nimport time\nfrom dataclasses import dataclass, field\nfrom typing import Protocol, TypeVar, Iterator, Callable\nfrom abc import ABC, abstractmethod\n\n\n# ============================================================\n# PART 1: Distribution Modeling Layer\n# ============================================================\n\n# --- Custom Exceptions ---\n\nclass DistributionError(Exception):\n """Base exception for distribution errors."""\n pass\n\n\nclass InvalidParameterError(DistributionError):\n """Raised when distribution parameters are invalid."""\n\n def __init__(self, param: str, value: float, reason: str) -> None:\n self.param = param\n self.value = value\n super().__init__(f"Invalid {param}={value}: {reason}")\n\n\n# --- Enums ---\n\nclass RiskLevel(enum.Enum):\n LOW = "low"\n MEDIUM = "medium"\n HIGH = "high"\n EXTREME = "extreme"\n\n @property\n def multiplier(self) -> float:\n multipliers = {\n RiskLevel.LOW: 0.5,\n RiskLevel.MEDIUM: 1.0,\n RiskLevel.HIGH: 1.5,\n RiskLevel.EXTREME: 2.0,\n }\n return multipliers[self]\n\n\n# --- Protocol (structural typing) ---\n\nT = TypeVar("T", bound="Distribution")\n\n\nclass Samplable(Protocol):\n def sample(self, n: int, *, seed: int | None = None) -> Iterator[float]: ...\n\n\n# --- Type Alias (3.12+ syntax) ---\n\ntype Samples = list[float]\ntype Parameters = dict[str, float]\n\n\n# --- Abstract Base Class ---\n\nclass Distribution(ABC):\n """Abstract base for probability distributions."""\n\n @abstractmethod\n def sample_one(self, rng: random.Random) -> float:\n """Generate a single sample."""\n ...\n\n @abstractmethod\n def expected_value(self) -> float:\n """Theoretical mean."""\n ...\n\n @abstractmethod\n def describe(self) -> str:\n """Human-readable description."""\n ...\n\n def sample(self, n: int, *, seed: int | None = None) -> Iterator[float]:\n """Generate n samples as a lazy iterator."""\n rng = random.Random(seed)\n for _ in itertools.repeat(None, n):\n yield self.sample_one(rng)\n\n def __repr__(self) -> str:\n return self.describe()\n\n\n# --- Frozen Dataclass Distributions ---\n\n@dataclass(frozen=True, slots=True)\nclass NormalReturn(Distribution):\n """Normally distributed returns."""\n\n mu: float\n sigma: float\n label: str = field(default="Normal", repr=False)\n\n def __post_init__(self) -> None:\n if self.sigma <= 0:\n raise InvalidParameterError("sigma", self.sigma, "must be positive")\n\n def sample_one(self, rng: random.Random) -> float:\n return rng.gauss(self.mu, self.sigma)\n\n def expected_value(self) -> float:\n return self.mu\n\n def describe(self) -> str:\n return f"{self.label}(\u03bc={self.mu:.2%}, \u03c3={self.sigma:.2%})"\n\n def __add__(self, other: NormalReturn) -> NormalReturn:\n """Sum of independent normals is normal."""\n combined_mu = self.mu + other.mu\n combined_sigma = math.sqrt(self.sigma**2 + other.sigma**2)\n return NormalReturn(combined_mu, combined_sigma, f"{self.label}+{other.label}")\n\n def __mul__(self, scalar: float) -> NormalReturn:\n """Scale a distribution."""\n return NormalReturn(self.mu * scalar, self.sigma * abs(scalar), self.label)\n\n def __rmul__(self, scalar: float) -> NormalReturn:\n return self.__mul__(scalar)\n\n\n@dataclass(frozen=True, slots=True)\nclass TriangularReturn(Distribution):\n """Triangular distribution for bounded uncertain returns."""\n\n low: float\n mode: float\n high: float\n label: str = field(default="Triangular", repr=False)\n\n def __post_init__(self) -> None:\n if not (self.low <= self.mode <= self.high):\n raise InvalidParameterError(\n "mode", self.mode,\n f"must satisfy low({self.low}) <= mode <= high({self.high})"\n )\n\n def sample_one(self, rng: random.Random) -> float:\n return rng.triangular(self.low, self.high, self.mode)\n\n def expected_value(self) -> float:\n return (self.low + self.mode + self.high) / 3\n\n def describe(self) -> str:\n return f"{self.label}(low={self.low:.2%}, mode={self.mode:.2%}, high={self.high:.2%})"\n\n\n@dataclass(frozen=True, slots=True)\nclass LogNormalReturn(Distribution):\n """Log-normal distribution for multiplicative returns."""\n\n mu: float\n sigma: float\n label: str = field(default="LogNormal", repr=False)\n\n def __post_init__(self) -> None:\n if self.sigma <= 0:\n raise InvalidParameterError("sigma", self.sigma, "must be positive")\n\n def sample_one(self, rng: random.Random) -> float:\n return math.exp(rng.gauss(self.mu, self.sigma)) - 1\n\n def expected_value(self) -> float:\n return math.exp(self.mu + self.sigma**2 / 2) - 1\n\n def describe(self) -> str:\n return f"{self.label}(\u03bc={self.mu:.3f}, \u03c3={self.sigma:.3f})"\n\n\n# --- @singledispatch factory ---\n\n@functools.singledispatch\ndef create_distribution(spec) -> Distribution:\n """Factory: create distribution from various input types."""\n raise TypeError(f"Cannot create distribution from {type(spec).__name__}")\n\n\n@create_distribution.register(dict)\ndef _from_dict(spec: dict) -> Distribution:\n match spec:\n case {"type": "normal", "mu": mu, "sigma": sigma, **rest}:\n return NormalReturn(mu, sigma, rest.get("label", "Normal"))\n case {"type": "triangular", "low": low, "mode": mode, "high": high, **rest}:\n return TriangularReturn(low, mode, high, rest.get("label", "Triangular"))\n case {"type": "lognormal", "mu": mu, "sigma": sigma, **rest}:\n return LogNormalReturn(mu, sigma, rest.get("label", "LogNormal"))\n case _:\n raise DistributionError(f"Unknown distribution spec: {spec}")\n\n\n@create_distribution.register(tuple)\ndef _from_tuple(spec: tuple) -> Distribution:\n """Quick creation: (low, mode, high) -> Triangular."""\n if len(spec) == 3:\n return TriangularReturn(*spec)\n raise DistributionError(f"Tuple must have 3 elements, got {len(spec)}")Your distribution layer is complete. You have a clean type hierarchy, immutable data models, operator overloading, and a flexible factory pattern. Next up, you'll build the simulation engine that uses these distributions to stream thousands of trial returns.
Build the Simulation Engine
Your distribution modeling layer gives you the building blocks for representing uncertain returns. Now you need an engine that can actually run simulations with those distributions.
In this step, you'll build the simulation engine layer. This includes a Strategy dataclass to pair distributions with investment metadata, decorators that add cross-cutting behavior, context managers for reproducible runs, and generator functions that stream simulation results lazily without loading thousands of trials into memory at once.
In this step, get ready to:
- Create a Strategy dataclass and three decorator patterns (simple, validation, and factory).
- Implement two styles of context manager for reproducible simulation control.
- Write generator-based simulation functions that stream results using itertools.
Create the Strategy dataclass and decorator patterns
Each investment strategy pairs a probability distribution with metadata like its name and risk level. The Strategy dataclass uses several advanced features to make strategies automatically sortable and keyword-only.
- In monte_carlo_engine.py, scroll to the bottom of the file (after the _from_tuple function).
- Add the PART 2 section header and the Strategy dataclass by copying this code below:
# ============================================================
# PART 2: Simulation Engine
# ============================================================
# --- Investment Strategy Dataclass ---
@dataclass(order=True, kw_only=True)
class Strategy:
"""An investment strategy with uncertain returns."""
sort_index: float = field(init=False, repr=False)
name: str
distribution: Distribution = field(compare=False)
risk_level: RiskLevel = field(compare=False)
allocation: float = field(default=0.0, compare=False)
def __post_init__(self) -> None:
self.sort_index = -self.distribution.expected_value()
@property
def allocation_pct(self) -> str:
return f"{self.allocation:.1%}"What does this code do?
This dataclass uses several advanced configuration options that control how instances behave.
- order=True auto-generates comparison dunders (__lt__, __le__, etc.) so strategies can be sorted.
- kw_only=True forces all fields to be passed as keyword arguments, preventing positional mix-ups.
- sort_index is computed in __post_init__ as the negative expected value, so sorting goes from highest to lowest return.
- field(compare=False) excludes fields from the auto-generated comparison, so only sort_index determines ordering.
Next, you need decorators to add cross-cutting behavior to functions. A simple decorator wraps a function with extra logic. A decorator factory takes arguments and returns a decorator.
- Add the timer and validation decorators below the Strategy class:
# --- Decorators ---
def timer(func: Callable) -> Callable:
"""Decorator to time function execution."""
@functools.wraps(func)
def wrapper(*args, **kwargs):
start = time.perf_counter()
result = func(*args, **kwargs)
elapsed = time.perf_counter() - start
print(f" \u23f1 {func.__name__} completed in {elapsed:.3f}s")
return result
return wrapper
def validate_allocation(func: Callable) -> Callable:
"""Decorator to validate allocations sum to 1.0."""
@functools.wraps(func)
def wrapper(allocations: dict[str, float], *args, **kwargs):
total = sum(allocations.values())
if not math.isclose(total, 1.0, rel_tol=1e-6):
raise ValueError(f"Allocations must sum to 1.0, got {total:.6f}")
return func(allocations, *args, **kwargs)
return wrapperWhat does this code do?
- @functools.wraps(func) preserves the original function's name and docstring on the wrapper. Without it, debugging becomes confusing because every decorated function would show up as "wrapper".
- timer measures execution time using time.perf_counter() (the highest-resolution clock) and prints the elapsed time.
- validate_allocation checks that portfolio weights sum to 1.0 before allowing the wrapped function to run. This is a guard pattern.
The third decorator is more complex. It's a decorator factory that accepts configuration arguments and returns the actual decorator.
- Add the retry decorator factory below validate_allocation:
def retry(max_attempts: int = 3, exceptions: tuple = (Exception,)):
"""Decorator factory with arguments for retry logic."""
def decorator(func: Callable) -> Callable:
@functools.wraps(func)
def wrapper(*args, **kwargs):
for attempt in range(1, max_attempts + 1):
try:
return func(*args, **kwargs)
except exceptions as e:
if attempt == max_attempts:
raise
print(f" \u26a0 Attempt {attempt} failed: {e}. Retrying...")
return None
return wrapper
return decoratorWhy three levels of nesting?
A decorator factory needs three layers: the outer function (retry) receives configuration, the middle function (decorator) receives the target function, and the inner function (wrapper) runs on each call.
This lets you write @retry(max_attempts=5) to customize behavior per function. The loop retries on specified exceptions and re-raises only after exhausting all attempts.
- Save monte_carlo_engine.py.
- Confirm there are no syntax errors by running:
python3 -c "import py_compile; py_compile.compile('monte_carlo_engine.py', doraise=True)"You should see no output, which means the file compiled without errors.
Seeing a SyntaxError?
- Check that each decorator function is properly closed with a matching return statement at the right indentation level.
- Make sure the retry function has three levels of def nesting, each indented one level deeper.
- Verify that f-strings use escaped braces correctly.
Help me debug a SyntaxError in my decorator code
Implement context managers for simulation control
Context managers give you setup/teardown guarantees using Python's with statement. You'll build two styles: a class-based context manager for complex state, and a generator-based one for lighter-weight tracking.
- Add the SimulationContext class below the retry decorator factory:
# --- Context Managers ---
class SimulationContext:
"""Context manager for reproducible simulation runs."""
def __init__(self, name: str, seed: int, num_trials: int) -> None:
self.name = name
self.seed = seed
self.num_trials = num_trials
self.rng: random.Random | None = None
self.results: list[float] = []
def __enter__(self) -> SimulationContext:
self.rng = random.Random(self.seed)
self.results = []
print(f"\n{'='*60}")
print(f" Simulation: {self.name}")
print(f" Trials: {self.num_trials:,} | Seed: {self.seed}")
print(f"{'='*60}")
return self
def __exit__(self, exc_type, exc_val, exc_tb) -> bool:
if exc_type is None:
print(f" \u2713 Completed {len(self.results):,} trials successfully")
else:
print(f" \u2717 Failed: {exc_val}")
return FalseWhat does this code do?
- __enter__ initializes a seeded RNG and prints a simulation banner. Returning self makes the context object available in the as clause.
- __exit__ runs cleanup logic. It checks whether an exception occurred and prints a success or failure message. Returning False means exceptions propagate normally.
- The fixed seed ensures every simulation run with the same seed produces identical results, which is critical for reproducibility.
The generator-based approach uses @contextlib.contextmanager to turn a generator function into a context manager with less boilerplate.
- Add the progress_tracker function below SimulationContext:
@contextlib.contextmanager
def progress_tracker(total: int, label: str = "Progress"):
"""Generator-based context manager that tracks simulation progress."""
state = {"completed": 0, "total": total}
def update(n: int = 1) -> None:
state["completed"] += n
if state["completed"] % max(total // 5, 1) == 0:
pct = state["completed"] / state["total"]
print(f" {label}: {pct:.0%}")
try:
yield update
finally:
print(f" {label}: Done ({state['completed']:,}/{state['total']:,})")What's happening here?
The @contextlib.contextmanager decorator turns a generator into a context manager. Everything before yield runs on entry, the yielded value becomes the as target, and the finally block runs on exit.
The yielded update function uses a mutable dict (state) as a closure variable to track progress. It prints every 20% of completion.
- Save monte_carlo_engine.py.
You won't see a visible change yet. The next chunk adds the simulation functions that use these context managers.
Write the simulation generator and streaming statistics
The simulation engine's core is a generator that streams portfolio returns trial-by-trial. Generators let you process results lazily, meaning you never need to hold all 10,000 trials in memory at once if you don't want to.
- Add the simulate_portfolio() generator below progress_tracker:
# --- Core Simulation Generator ---
def simulate_portfolio(
strategies: list[Strategy],
allocations: dict[str, float],
num_trials: int,
*,
seed: int = 42,
) -> Iterator[float]:
"""Generator yielding portfolio returns for each trial."""
rng = random.Random(seed)
alloc_pairs = [(s, allocations[s.name]) for s in strategies]
for _ in itertools.repeat(None, num_trials):
trial_return = sum(
alloc * s.distribution.sample_one(rng)
for s, alloc in alloc_pairs
)
yield trial_returnWhat does this code do?
- itertools.repeat(None, num_trials) creates an efficient loop counter without allocating a list. It's faster than range() for large iteration counts where you don't need the index.
- The * in the function signature makes seed keyword-only, preventing accidental positional passing.
- Each yield pauses the function and produces one trial's weighted portfolio return. The caller controls how many results to consume.
The final piece computes running statistics from a stream of returns using itertools composition.
- Add compute_running_stats() below simulate_portfolio():
# --- Running Statistics via itertools ---
def compute_running_stats(returns: Iterator[float]) -> Iterator[tuple]:
"""Compute running statistics using itertools patterns."""
t0, t1, t2 = itertools.tee(returns, 3)
running_mean = map(
lambda pair: pair[0] / pair[1],
zip(itertools.accumulate(t0), itertools.count(1))
)
running_max_iter = itertools.accumulate(t1, func=max)
running_min_iter = itertools.accumulate(t2, func=min)
return zip(itertools.count(1), running_mean, running_min_iter, running_max_iter)What does this code do?
- itertools.tee(returns, 3) creates three independent iterators from one input stream. Each can be consumed independently without affecting the others.
- itertools.accumulate(t0) produces a running sum. Dividing by itertools.count(1) (1, 2, 3...) gives you the running mean at each step.
- accumulate(t1, func=max) tracks the running maximum by always keeping the larger of the accumulated value and the new element.
- The final zip combines all streams into tuples of (trial_number, running_mean, running_min, running_max). Everything stays lazy until consumed.
- Save monte_carlo_engine.py.
- Verify the full file has no errors by running:
python3 -c "import py_compile; py_compile.compile('monte_carlo_engine.py', doraise=True)"No output means the entire simulation engine layer compiles cleanly.
Seeing an IndentationError or NameError?
- Check that compute_running_stats is at the top-level indentation (no extra indent from a class or function above it).
- Confirm you have import itertools and import contextlib at the top of the file from Step 1.
- Make sure the lambda inside running_mean has the colon and expression on the same line.
Help me fix an error in my simulation engine code
✔️ Awesome, I've got everything!
Great. Double check you've saved your file before moving on.
ⓧ I'd like to double check the full code
Here's the full monte_carlo_engine.py as it should look at the end of this step:
"""
Monte Carlo Decision Analysis Engine
Pure Python Standard Library Implementation
"""
from __future__ import annotations
import random
import statistics
import itertools
import functools
import contextlib
import collections
import enum
import math
import time
from dataclasses import dataclass, field
from typing import Protocol, TypeVar, Iterator, Callable
from abc import ABC, abstractmethod
# ============================================================
# PART 1: Distribution Modeling Layer
# ============================================================
class DistributionError(Exception):
pass
class InvalidParameterError(DistributionError):
def __init__(self, param: str, value: float, reason: str) -> None:
self.param = param
self.value = value
super().__init__(f"Invalid {param}={value}: {reason}")
class RiskLevel(enum.Enum):
LOW = "low"
MEDIUM = "medium"
HIGH = "high"
EXTREME = "extreme"
@property
def multiplier(self) -> float:
multipliers = {
RiskLevel.LOW: 0.5,
RiskLevel.MEDIUM: 1.0,
RiskLevel.HIGH: 1.5,
RiskLevel.EXTREME: 2.0,
}
return multipliers[self]
T = TypeVar("T", bound="Distribution")
class Samplable(Protocol):
def sample(self, n: int, *, seed: int | None = None) -> Iterator[float]: ...
type Samples = list[float]
type Parameters = dict[str, float]
class Distribution(ABC):
@abstractmethod
def sample_one(self, rng: random.Random) -> float: ...
@abstractmethod
def expected_value(self) -> float: ...
@abstractmethod
def describe(self) -> str: ...
def sample(self, n: int, *, seed: int | None = None) -> Iterator[float]:
rng = random.Random(seed)
for _ in itertools.repeat(None, n):
yield self.sample_one(rng)
def __repr__(self) -> str:
return self.describe()
@dataclass(frozen=True, slots=True)
class NormalReturn(Distribution):
mu: float
sigma: float
label: str = field(default="Normal", repr=False)
def __post_init__(self) -> None:
if self.sigma <= 0:
raise InvalidParameterError("sigma", self.sigma, "must be positive")
def sample_one(self, rng: random.Random) -> float:
return rng.gauss(self.mu, self.sigma)
def expected_value(self) -> float:
return self.mu
def describe(self) -> str:
return f"{self.label}(\u03bc={self.mu:.2%}, \u03c3={self.sigma:.2%})"
def __add__(self, other: NormalReturn) -> NormalReturn:
combined_mu = self.mu + other.mu
combined_sigma = math.sqrt(self.sigma**2 + other.sigma**2)
return NormalReturn(combined_mu, combined_sigma, f"{self.label}+{other.label}")
def __mul__(self, scalar: float) -> NormalReturn:
return NormalReturn(self.mu * scalar, self.sigma * abs(scalar), self.label)
def __rmul__(self, scalar: float) -> NormalReturn:
return self.__mul__(scalar)
@dataclass(frozen=True, slots=True)
class TriangularReturn(Distribution):
low: float
mode: float
high: float
label: str = field(default="Triangular", repr=False)
def __post_init__(self) -> None:
if not (self.low <= self.mode <= self.high):
raise InvalidParameterError("mode", self.mode, f"must satisfy low({self.low}) <= mode <= high({self.high})")
def sample_one(self, rng: random.Random) -> float:
return rng.triangular(self.low, self.high, self.mode)
def expected_value(self) -> float:
return (self.low + self.mode + self.high) / 3
def describe(self) -> str:
return f"{self.label}(low={self.low:.2%}, mode={self.mode:.2%}, high={self.high:.2%})"
@dataclass(frozen=True, slots=True)
class LogNormalReturn(Distribution):
mu: float
sigma: float
label: str = field(default="LogNormal", repr=False)
def __post_init__(self) -> None:
if self.sigma <= 0:
raise InvalidParameterError("sigma", self.sigma, "must be positive")
def sample_one(self, rng: random.Random) -> float:
return math.exp(rng.gauss(self.mu, self.sigma)) - 1
def expected_value(self) -> float:
return math.exp(self.mu + self.sigma**2 / 2) - 1
def describe(self) -> str:
return f"{self.label}(\u03bc={self.mu:.3f}, \u03c3={self.sigma:.3f})"
@functools.singledispatch
def create_distribution(spec) -> Distribution:
raise TypeError(f"Cannot create distribution from {type(spec).__name__}")
@create_distribution.register(dict)
def _from_dict(spec: dict) -> Distribution:
match spec:
case {"type": "normal", "mu": mu, "sigma": sigma, **rest}:
return NormalReturn(mu, sigma, rest.get("label", "Normal"))
case {"type": "triangular", "low": low, "mode": mode, "high": high, **rest}:
return TriangularReturn(low, mode, high, rest.get("label", "Triangular"))
case {"type": "lognormal", "mu": mu, "sigma": sigma, **rest}:
return LogNormalReturn(mu, sigma, rest.get("label", "LogNormal"))
case _:
raise DistributionError(f"Unknown distribution spec: {spec}")
@create_distribution.register(tuple)
def _from_tuple(spec: tuple) -> Distribution:
if len(spec) == 3:
return TriangularReturn(*spec)
raise DistributionError(f"Tuple must have 3 elements, got {len(spec)}")
# ============================================================
# PART 2: Simulation Engine
# ============================================================
@dataclass(order=True, kw_only=True)
class Strategy:
sort_index: float = field(init=False, repr=False)
name: str
distribution: Distribution = field(compare=False)
risk_level: RiskLevel = field(compare=False)
allocation: float = field(default=0.0, compare=False)
def __post_init__(self) -> None:
self.sort_index = -self.distribution.expected_value()
@property
def allocation_pct(self) -> str:
return f"{self.allocation:.1%}"
def timer(func: Callable) -> Callable:
@functools.wraps(func)
def wrapper(*args, **kwargs):
start = time.perf_counter()
result = func(*args, **kwargs)
elapsed = time.perf_counter() - start
print(f" \u23f1 {func.__name__} completed in {elapsed:.3f}s")
return result
return wrapper
def validate_allocation(func: Callable) -> Callable:
@functools.wraps(func)
def wrapper(allocations: dict[str, float], *args, **kwargs):
total = sum(allocations.values())
if not math.isclose(total, 1.0, rel_tol=1e-6):
raise ValueError(f"Allocations must sum to 1.0, got {total:.6f}")
return func(allocations, *args, **kwargs)
return wrapper
def retry(max_attempts: int = 3, exceptions: tuple = (Exception,)):
def decorator(func: Callable) -> Callable:
@functools.wraps(func)
def wrapper(*args, **kwargs):
for attempt in range(1, max_attempts + 1):
try:
return func(*args, **kwargs)
except exceptions as e:
if attempt == max_attempts:
raise
print(f" \u26a0 Attempt {attempt} failed: {e}. Retrying...")
return None
return wrapper
return decorator
class SimulationContext:
def __init__(self, name: str, seed: int, num_trials: int) -> None:
self.name = name
self.seed = seed
self.num_trials = num_trials
self.rng: random.Random | None = None
self.results: list[float] = []
def __enter__(self) -> SimulationContext:
self.rng = random.Random(self.seed)
self.results = []
print(f"\n{'='*60}")
print(f" Simulation: {self.name}")
print(f" Trials: {self.num_trials:,} | Seed: {self.seed}")
print(f"{'='*60}")
return self
def __exit__(self, exc_type, exc_val, exc_tb) -> bool:
if exc_type is None:
print(f" \u2713 Completed {len(self.results):,} trials successfully")
else:
print(f" \u2717 Failed: {exc_val}")
return False
@contextlib.contextmanager
def progress_tracker(total: int, label: str = "Progress"):
state = {"completed": 0, "total": total}
def update(n: int = 1) -> None:
state["completed"] += n
if state["completed"] % max(total // 5, 1) == 0:
pct = state["completed"] / state["total"]
print(f" {label}: {pct:.0%}")
try:
yield update
finally:
print(f" {label}: Done ({state['completed']:,}/{state['total']:,})")
def simulate_portfolio(
strategies: list[Strategy],
allocations: dict[str, float],
num_trials: int,
*,
seed: int = 42,
) -> Iterator[float]:
rng = random.Random(seed)
alloc_pairs = [(s, allocations[s.name]) for s in strategies]
for _ in itertools.repeat(None, num_trials):
trial_return = sum(
alloc * s.distribution.sample_one(rng)
for s, alloc in alloc_pairs
)
yield trial_return
def compute_running_stats(returns: Iterator[float]) -> Iterator[tuple]:
t0, t1, t2 = itertools.tee(returns, 3)
running_mean = map(
lambda pair: pair[0] / pair[1],
zip(itertools.accumulate(t0), itertools.count(1))
)
running_max_iter = itertools.accumulate(t1, func=max)
running_min_iter = itertools.accumulate(t2, func=min)
return zip(itertools.count(1), running_mean, running_min_iter, running_max_iter)Your simulation engine is ready to run trials. Next up, you'll build the decision analysis framework that uses this engine to search through thousands of allocation combinations and find the optimal investment split.
Build the Decision Analysis Framework
Your simulation engine can stream thousands of portfolio returns. Now you need a way to make sense of those results and find the best possible allocation.
In this step, you will build a decision analysis framework that evaluates portfolio performance using cached properties, searches through every valid allocation combination using combinatoric generation, and selects an optimization objective using structural pattern matching.
In this step, get ready to:
- Create a PortfolioResult dataclass with lazy-computed risk metrics.
- Build an allocation generator that enumerates every valid portfolio split.
- Implement the optimizer that finds the best allocation using pattern matching.
Create the PortfolioResult dataclass with cached metrics
Each simulation run produces a list of returns. You need a structured way to compute and cache statistics from those results without recalculating them on every access.
@functools.cached_property computes a value once on first access, then stores it as a regular attribute. This is ideal for expensive statistical computations you reference multiple times.
- In monte_carlo_engine.py, scroll to the bottom of your file (after the compute_running_stats function).
- Add the Part 3 section header and the first half of the PortfolioResult class by copying the following code:
# ============================================================
# PART 3: Decision Analysis Framework
# ============================================================
@dataclass
class PortfolioResult:
allocation: dict[str, float]
returns: Samples
_stats_cache: dict = field(default_factory=dict, repr=False)
@functools.cached_property
def mean_return(self) -> float:
return statistics.fmean(self.returns)
@functools.cached_property
def volatility(self) -> float:
return statistics.stdev(self.returns)
@functools.cached_property
def sharpe_ratio(self) -> float:
# Walrus operator assigns volatility and checks for zero in one line
if (vol := self.volatility) == 0:
return 0.0
return self.mean_return / volWhat does this code do?
@functools.cached_property makes each metric compute only once. The first time you access result.mean_return, Python calls the method. Every subsequent access returns the cached value instantly.
The walrus operator := in sharpe_ratio assigns self.volatility to vol and checks if it equals zero, all in one expression. This avoids calling self.volatility twice.
The class needs risk metrics and comparison methods too. These let you rank portfolios and identify tail risk.
- Add the remaining cached properties and dunder methods directly below the sharpe_ratio property, still inside the PortfolioResult class:
@functools.cached_property
def var_95(self) -> float:
# Value at Risk: the 5th percentile of returns
sorted_returns = sorted(self.returns)
idx = int(len(sorted_returns) * 0.05)
return sorted_returns[idx]
@functools.cached_property
def cvar_95(self) -> float:
# Conditional VaR: average of returns below the VaR threshold
sorted_returns = sorted(self.returns)
idx = int(len(sorted_returns) * 0.05)
tail = sorted_returns[:idx] if idx > 0 else sorted_returns[:1]
return statistics.fmean(tail)
@functools.cached_property
def quantiles(self) -> list[float]:
return statistics.quantiles(self.returns, n=4)
def probability_of_loss(self) -> float:
# Fraction of trials that produced negative returns
return sum(1 for r in self.returns if r < 0) / len(self.returns)
def __lt__(self, other: PortfolioResult) -> bool:
return self.sharpe_ratio < other.sharpe_ratio
def __eq__(self, other: object) -> bool:
if not isinstance(other, PortfolioResult):
return NotImplemented
return math.isclose(self.sharpe_ratio, other.sharpe_ratio, rel_tol=1e-6)What do these risk metrics mean?
- var_95 (Value at Risk) answers: "In the worst 5% of scenarios, how much do I lose at minimum?"
- cvar_95 (Conditional VaR) answers: "When things go bad, how bad is the average loss?"
- __lt__ and __eq__ let Python compare PortfolioResult objects using < and == operators, comparing by Sharpe ratio.
- Save monte_carlo_engine.py.
- You won't see a change yet. The next chunk uses this class to evaluate allocations.
Seeing an indentation error?
Every method from mean_return through __eq__ must be indented one level inside the PortfolioResult class. Make sure you pasted the second block at the same indent level as sharpe_ratio, not outside the class.
help me fix indentation in my PortfolioResult class
Build the allocation generator
To find the optimal allocation, you need to try every valid combination of weights across three strategies. If you allocate in 10% increments, that means every combination of 0%, 10%, 20%... up to 100% that sums to exactly 1.0.
- Add the generate_allocations() generator function below the PortfolioResult class (outside the class, at the module level):
def generate_allocations(
strategy_names: list[str],
step: float = 0.1,
) -> Iterator[dict[str, float]]:
# Number of possible weight values (0.0, 0.1, 0.2, ... 1.0)
n = len(strategy_names)
steps = int(1.0 / step) + 1
possible_values = [round(i * step, 2) for i in range(steps)]
# itertools.product generates every combination of n weights
for combo in itertools.product(possible_values, repeat=n):
# Only yield combinations that sum to 1.0
if math.isclose(sum(combo), 1.0, rel_tol=1e-9):
yield dict(zip(strategy_names, combo))What does this code do?
- itertools.product(possible_values, repeat=n) generates the Cartesian product of all weight values for each strategy. With 3 strategies and 11 possible values (0.0 to 1.0 in 0.1 steps), that is 11³ = 1,331 combinations to test.
- math.isclose(sum(combo), 1.0) filters to only valid allocations where all weights sum to exactly 100%. This yields 66 valid combinations for 3 strategies.
- dict(zip(strategy_names, combo)) pairs each strategy name with its weight, yielding a dictionary like {"Early-Stage": 0.3, "Growth": 0.4, "RealEstate": 0.3}.
- Save monte_carlo_engine.py.
- You won't see output yet. The optimizer in the next chunk consumes this generator.
Getting a NameError on Iterator?
Make sure your imports at the top of the file include from typing import Protocol, TypeVar, Iterator, Callable. The Iterator type hint was imported in Step 1.
help me fix NameError in my generate_allocations function
Implement the optimization function with pattern matching
The optimizer needs to support multiple objectives. A fund manager might optimize for the highest Sharpe ratio, the best raw return, the lowest volatility, or the best Value at Risk. Structural pattern matching makes this selection clean and extensible.
- Add the find_optimal_allocation() function below generate_allocations():
@timer
def find_optimal_allocation(
strategies: list[Strategy],
num_trials: int = 10_000,
step: float = 0.1,
objective: str = "sharpe",
*,
seed: int = 42,
) -> PortfolioResult:
strategy_names = [s.name for s in strategies]
# Structural pattern matching selects the scoring function
match objective:
case "sharpe":
key_func = lambda r: r.sharpe_ratio
case "return":
key_func = lambda r: r.mean_return
case "risk":
key_func = lambda r: -r.volatility
case "var":
key_func = lambda r: r.var_95
case str() as unknown:
raise ValueError(f"Unknown objective: {unknown}")
best: PortfolioResult | None = None
tested = 0
for allocation in generate_allocations(strategy_names, step):
returns = list(simulate_portfolio(strategies, allocation, num_trials, seed=seed))
result = PortfolioResult(allocation=allocation, returns=returns)
# Walrus operator assigns score and compares in one expression
if best is None or (score := key_func(result)) > key_func(best):
best = result
tested += 1
print(f" Tested {tested} improving allocations")
return bestWhat does this code do?
- The match objective block uses structural pattern matching to select a scoring function based on the objective string. The case str() as unknown pattern catches any unrecognized string and binds it to unknown for the error message.
- The @timer decorator (built in Step 3) wraps this function to print its execution time. With 66 allocations and 10,000 trials each, this runs 660,000 simulations.
- The walrus operator (score := key_func(result)) > key_func(best) computes the score, saves it to score, and compares it against the current best, all in one line.
Why pattern matching instead of if/elif?
Pattern matching provides exhaustiveness checking and structured decomposition that if/elif chains cannot. The case str() as unknown pattern matches any string type while binding its value. This makes the catch-all case both type-safe and informative.
- Save monte_carlo_engine.py.
- You won't see output yet. The next step connects this framework to the reporting layer and runs the full engine.
Seeing a SyntaxError on match/case?
Structural pattern matching requires Python 3.10 or higher. Run python3 --version to confirm your version. If you see 3.9 or lower, you need to upgrade to 3.12+ as required by this project.
help me fix a SyntaxError with match/case in my Python file
✔️ Awesome, I've got everything!
Great. Double check you have saved your file before moving on.
ⓧ I'd like to double check the full code
Here is the full monte_carlo_engine.py as it should look at the end of this step:
"""
Monte Carlo Decision Analysis Engine
Pure Python Standard Library Implementation
"""
from __future__ import annotations
import random
import statistics
import itertools
import functools
import contextlib
import collections
import enum
import math
import time
from dataclasses import dataclass, field
from typing import Protocol, TypeVar, Iterator, Callable
from abc import ABC, abstractmethod
# ============================================================
# PART 1: Distribution Modeling Layer
# ============================================================
class DistributionError(Exception):
pass
class InvalidParameterError(DistributionError):
def __init__(self, param: str, value: float, reason: str) -> None:
self.param = param
self.value = value
super().__init__(f"Invalid {param}={value}: {reason}")
class RiskLevel(enum.Enum):
LOW = "low"
MEDIUM = "medium"
HIGH = "high"
EXTREME = "extreme"
@property
def multiplier(self) -> float:
multipliers = {
RiskLevel.LOW: 0.5,
RiskLevel.MEDIUM: 1.0,
RiskLevel.HIGH: 1.5,
RiskLevel.EXTREME: 2.0,
}
return multipliers[self]
T = TypeVar("T", bound="Distribution")
class Samplable(Protocol):
def sample(self, n: int, *, seed: int | None = None) -> Iterator[float]: ...
type Samples = list[float]
type Parameters = dict[str, float]
class Distribution(ABC):
@abstractmethod
def sample_one(self, rng: random.Random) -> float: ...
@abstractmethod
def expected_value(self) -> float: ...
@abstractmethod
def describe(self) -> str: ...
def sample(self, n: int, *, seed: int | None = None) -> Iterator[float]:
rng = random.Random(seed)
for _ in itertools.repeat(None, n):
yield self.sample_one(rng)
def __repr__(self) -> str:
return self.describe()
@dataclass(frozen=True, slots=True)
class NormalReturn(Distribution):
mu: float
sigma: float
label: str = field(default="Normal", repr=False)
def __post_init__(self) -> None:
if self.sigma <= 0:
raise InvalidParameterError("sigma", self.sigma, "must be positive")
def sample_one(self, rng: random.Random) -> float:
return rng.gauss(self.mu, self.sigma)
def expected_value(self) -> float:
return self.mu
def describe(self) -> str:
return f"{self.label}(\u03bc={self.mu:.2%}, \u03c3={self.sigma:.2%})"
def __add__(self, other: NormalReturn) -> NormalReturn:
combined_mu = self.mu + other.mu
combined_sigma = math.sqrt(self.sigma**2 + other.sigma**2)
return NormalReturn(combined_mu, combined_sigma, f"{self.label}+{other.label}")
def __mul__(self, scalar: float) -> NormalReturn:
return NormalReturn(self.mu * scalar, self.sigma * abs(scalar), self.label)
def __rmul__(self, scalar: float) -> NormalReturn:
return self.__mul__(scalar)
@dataclass(frozen=True, slots=True)
class TriangularReturn(Distribution):
low: float
mode: float
high: float
label: str = field(default="Triangular", repr=False)
def __post_init__(self) -> None:
if not (self.low <= self.mode <= self.high):
raise InvalidParameterError("mode", self.mode, f"must satisfy low({self.low}) <= mode <= high({self.high})")
def sample_one(self, rng: random.Random) -> float:
return rng.triangular(self.low, self.high, self.mode)
def expected_value(self) -> float:
return (self.low + self.mode + self.high) / 3
def describe(self) -> str:
return f"{self.label}(low={self.low:.2%}, mode={self.mode:.2%}, high={self.high:.2%})"
@dataclass(frozen=True, slots=True)
class LogNormalReturn(Distribution):
mu: float
sigma: float
label: str = field(default="LogNormal", repr=False)
def __post_init__(self) -> None:
if self.sigma <= 0:
raise InvalidParameterError("sigma", self.sigma, "must be positive")
def sample_one(self, rng: random.Random) -> float:
return math.exp(rng.gauss(self.mu, self.sigma)) - 1
def expected_value(self) -> float:
return math.exp(self.mu + self.sigma**2 / 2) - 1
def describe(self) -> str:
return f"{self.label}(\u03bc={self.mu:.3f}, \u03c3={self.sigma:.3f})"
@functools.singledispatch
def create_distribution(spec) -> Distribution:
raise TypeError(f"Cannot create distribution from {type(spec).__name__}")
@create_distribution.register(dict)
def _from_dict(spec: dict) -> Distribution:
match spec:
case {"type": "normal", "mu": mu, "sigma": sigma, **rest}:
return NormalReturn(mu, sigma, rest.get("label", "Normal"))
case {"type": "triangular", "low": low, "mode": mode, "high": high, **rest}:
return TriangularReturn(low, mode, high, rest.get("label", "Triangular"))
case {"type": "lognormal", "mu": mu, "sigma": sigma, **rest}:
return LogNormalReturn(mu, sigma, rest.get("label", "LogNormal"))
case _:
raise DistributionError(f"Unknown distribution spec: {spec}")
@create_distribution.register(tuple)
def _from_tuple(spec: tuple) -> Distribution:
if len(spec) == 3:
return TriangularReturn(*spec)
raise DistributionError(f"Tuple must have 3 elements, got {len(spec)}")
# ============================================================
# PART 2: Simulation Engine
# ============================================================
@dataclass(order=True, kw_only=True)
class Strategy:
sort_index: float = field(init=False, repr=False)
name: str
distribution: Distribution = field(compare=False)
risk_level: RiskLevel = field(compare=False)
allocation: float = field(default=0.0, compare=False)
def __post_init__(self) -> None:
self.sort_index = -self.distribution.expected_value()
@property
def allocation_pct(self) -> str:
return f"{self.allocation:.1%}"
def timer(func: Callable) -> Callable:
@functools.wraps(func)
def wrapper(*args, **kwargs):
start = time.perf_counter()
result = func(*args, **kwargs)
elapsed = time.perf_counter() - start
print(f" \u23f1 {func.__name__} completed in {elapsed:.3f}s")
return result
return wrapper
def validate_allocation(func: Callable) -> Callable:
@functools.wraps(func)
def wrapper(allocations: dict[str, float], *args, **kwargs):
total = sum(allocations.values())
if not math.isclose(total, 1.0, rel_tol=1e-6):
raise ValueError(f"Allocations must sum to 1.0, got {total:.6f}")
return func(allocations, *args, **kwargs)
return wrapper
def retry(max_attempts: int = 3, exceptions: tuple = (Exception,)):
def decorator(func: Callable) -> Callable:
@functools.wraps(func)
def wrapper(*args, **kwargs):
for attempt in range(1, max_attempts + 1):
try:
return func(*args, **kwargs)
except exceptions as e:
if attempt == max_attempts:
raise
print(f" \u26a0 Attempt {attempt} failed: {e}. Retrying...")
return None
return wrapper
return decorator
class SimulationContext:
def __init__(self, name: str, seed: int, num_trials: int) -> None:
self.name = name
self.seed = seed
self.num_trials = num_trials
self.rng: random.Random | None = None
self.results: list[float] = []
def __enter__(self) -> SimulationContext:
self.rng = random.Random(self.seed)
self.results = []
print(f"\n{'='*60}")
print(f" Simulation: {self.name}")
print(f" Trials: {self.num_trials:,} | Seed: {self.seed}")
print(f"{'='*60}")
return self
def __exit__(self, exc_type, exc_val, exc_tb) -> bool:
if exc_type is None:
print(f" \u2713 Completed {len(self.results):,} trials successfully")
else:
print(f" \u2717 Failed: {exc_val}")
return False
@contextlib.contextmanager
def progress_tracker(total: int, label: str = "Progress"):
state = {"completed": 0, "total": total}
def update(n: int = 1) -> None:
state["completed"] += n
if state["completed"] % max(total // 5, 1) == 0:
pct = state["completed"] / state["total"]
print(f" {label}: {pct:.0%}")
try:
yield update
finally:
print(f" {label}: Done ({state['completed']:,}/{state['total']:,})")
def simulate_portfolio(
strategies: list[Strategy],
allocations: dict[str, float],
num_trials: int,
*,
seed: int = 42,
) -> Iterator[float]:
rng = random.Random(seed)
alloc_pairs = [(s, allocations[s.name]) for s in strategies]
for _ in itertools.repeat(None, num_trials):
trial_return = sum(
alloc * s.distribution.sample_one(rng)
for s, alloc in alloc_pairs
)
yield trial_return
def compute_running_stats(returns: Iterator[float]) -> Iterator[tuple]:
t0, t1, t2 = itertools.tee(returns, 3)
running_mean = map(
lambda pair: pair[0] / pair[1],
zip(itertools.accumulate(t0), itertools.count(1))
)
running_max_iter = itertools.accumulate(t1, func=max)
running_min_iter = itertools.accumulate(t2, func=min)
return zip(itertools.count(1), running_mean, running_min_iter, running_max_iter)
# ============================================================
# PART 3: Decision Analysis Framework
# ============================================================
@dataclass
class PortfolioResult:
allocation: dict[str, float]
returns: Samples
_stats_cache: dict = field(default_factory=dict, repr=False)
@functools.cached_property
def mean_return(self) -> float:
return statistics.fmean(self.returns)
@functools.cached_property
def volatility(self) -> float:
return statistics.stdev(self.returns)
@functools.cached_property
def sharpe_ratio(self) -> float:
if (vol := self.volatility) == 0:
return 0.0
return self.mean_return / vol
@functools.cached_property
def var_95(self) -> float:
sorted_returns = sorted(self.returns)
idx = int(len(sorted_returns) * 0.05)
return sorted_returns[idx]
@functools.cached_property
def cvar_95(self) -> float:
sorted_returns = sorted(self.returns)
idx = int(len(sorted_returns) * 0.05)
tail = sorted_returns[:idx] if idx > 0 else sorted_returns[:1]
return statistics.fmean(tail)
@functools.cached_property
def quantiles(self) -> list[float]:
return statistics.quantiles(self.returns, n=4)
def probability_of_loss(self) -> float:
return sum(1 for r in self.returns if r < 0) / len(self.returns)
def __lt__(self, other: PortfolioResult) -> bool:
return self.sharpe_ratio < other.sharpe_ratio
def __eq__(self, other: object) -> bool:
if not isinstance(other, PortfolioResult):
return NotImplemented
return math.isclose(self.sharpe_ratio, other.sharpe_ratio, rel_tol=1e-6)
def generate_allocations(
strategy_names: list[str],
step: float = 0.1,
) -> Iterator[dict[str, float]]:
n = len(strategy_names)
steps = int(1.0 / step) + 1
possible_values = [round(i * step, 2) for i in range(steps)]
for combo in itertools.product(possible_values, repeat=n):
if math.isclose(sum(combo), 1.0, rel_tol=1e-9):
yield dict(zip(strategy_names, combo))
@timer
def find_optimal_allocation(
strategies: list[Strategy],
num_trials: int = 10_000,
step: float = 0.1,
objective: str = "sharpe",
*,
seed: int = 42,
) -> PortfolioResult:
strategy_names = [s.name for s in strategies]
match objective:
case "sharpe":
key_func = lambda r: r.sharpe_ratio
case "return":
key_func = lambda r: r.mean_return
case "risk":
key_func = lambda r: -r.volatility
case "var":
key_func = lambda r: r.var_95
case str() as unknown:
raise ValueError(f"Unknown objective: {unknown}")
best: PortfolioResult | None = None
tested = 0
for allocation in generate_allocations(strategy_names, step):
returns = list(simulate_portfolio(strategies, allocation, num_trials, seed=seed))
result = PortfolioResult(allocation=allocation, returns=returns)
if best is None or (score := key_func(result)) > key_func(best):
best = result
tested += 1
print(f" Tested {tested} improving allocations")
return bestYour decision analysis framework is complete. It can evaluate any portfolio allocation, cache expensive computations, and search through all valid combinations to find the best one. Next up, you will build the reporting layer and run the full engine end-to-end.
Build the Reporting System and Run the Engine
Your decision analysis framework can now find the optimal allocation across thousands of simulated portfolios. The next challenge is making those results readable and actionable.
In this step, you will build a reporting system that presents the optimal portfolio with statistical depth. You will also write the main() function that ties the entire engine together and runs it end-to-end.
In this step, get ready to:
- Build correlation matrix and report printing functions using the statistics module.
- Implement an ASCII histogram and a main function that orchestrates the full analysis.
- Run the complete engine and interpret its output.
Build the correlation matrix and report functions
A correlation matrix reveals how strategy returns move together. If two strategies are highly correlated, diversifying between them adds less protection than expected.
- In monte_carlo_engine.py, scroll to the bottom of your file (after the find_optimal_allocation() function).
- Add the Part 4 section header and the correlation matrix function by copying the following code:
# ============================================================
# PART 4: Reporting & Analysis
# ============================================================
def correlation_matrix(
strategies: list[Strategy],
n_samples: int = 5000,
seed: int = 42,
) -> dict[tuple[str, str], float]:
samples: dict[str, Samples] = {
s.name: list(s.distribution.sample(n_samples, seed=seed))
for s in strategies
}
correlations: dict[tuple[str, str], float] = {}
for (name_a, data_a), (name_b, data_b) in itertools.combinations(
samples.items(), 2
):
corr = statistics.correlation(data_a, data_b)
correlations[(name_a, name_b)] = corr
correlations[(name_b, name_a)] = corr
for name in samples:
correlations[(name, name)] = 1.0
return correlationsWhat does this code do?
This function generates sample data from each strategy's distribution, then uses itertools.combinations() to pair up every combination of strategies.
For each pair, statistics.correlation() computes the Pearson correlation coefficient. The function stores both directions (A,B) and (B,A) for easy lookup, plus 1.0 for self-correlation.
The report function produces the full human-readable output. It uses several powerful features from the statistics module to fit a normal distribution and compute tail probabilities.
- Add the print_report() function directly below correlation_matrix():
def print_report(
result: PortfolioResult,
strategies: list[Strategy],
fund_size: float = 10_000_000,
) -> None:
print(f"\n{'\u2550'*60}")
print(f" \U0001f4ca OPTIMAL PORTFOLIO ANALYSIS")
print(f"{'\u2550'*60}")
print(f"\n \U0001f4b0 Allocation (Fund: ${fund_size:,.0f}):")
for name, weight in result.allocation.items():
dollar_amount = weight * fund_size
print(f" {name:.<30} {weight:>6.1%} (${dollar_amount:>12,.0f})")
print(f"\n \U0001f4c8 Performance Metrics:")
print(f" Expected Return: {result.mean_return:>8.2%}")
print(f" Volatility (\u03c3): {result.volatility:>8.2%}")
print(f" Sharpe Ratio: {result.sharpe_ratio:>8.3f}")
print(f" Value at Risk (95%): {result.var_95:>8.2%}")
print(f" CVaR (95%): {result.cvar_95:>8.2%}")
print(f" P(Loss): {result.probability_of_loss():>8.2%}")
dist = statistics.NormalDist.from_samples(result.returns)
print(f"\n \U0001f4c9 Return Distribution (fitted Normal):")
print(f" Mean: {dist.mean:>8.2%}")
print(f" Stdev: {dist.stdev:>8.2%}")
print(f" P(return > 20%): {1 - dist.cdf(0.20):>8.2%}")
print(f" P(return < -10%): {dist.cdf(-0.10):>8.2%}")
q1, q2, q3 = result.quantiles
print(f"\n \U0001f4ca Quartiles:")
print(f" Q1 (25th): {q1:>8.2%}")
print(f" Q2 (50th): {q2:>8.2%}")
print(f" Q3 (75th): {q3:>8.2%}")
print(f"\n \U0001f4ca Return Distribution Histogram:")
_print_histogram(result.returns)
correlations = correlation_matrix(strategies)
names = [s.name for s in strategies]
print(f"\n \U0001f517 Correlation Matrix:")
header = " " + "".join(f"{n[:10]:>12}" for n in names)
print(header)
for name_a in names:
row = f" {name_a[:10]:<10}"
for name_b in names:
corr = correlations[(name_a, name_b)]
row += f"{corr:>12.3f}"
print(row)
print(f"\n{'\u2550'*60}")What does this report function do?
- statistics.NormalDist.from_samples() fits a normal distribution to the simulated returns, giving you access to .mean, .stdev, and .cdf() for tail probability calculations.
- F-string format specs like :.2% (percentage with 2 decimals), :>8 (right-align in 8 chars), and :,.0f (comma-separated integer) produce aligned, readable output.
- The correlation matrix section reuses the correlation_matrix() function you just wrote.
- Save monte_carlo_engine.py.
Seeing a syntax error in the f-strings?
The nested quotes inside f-strings (like f"{'═'*60}") require Python 3.12+. If you see a SyntaxError, confirm your Python version with python3 --version.
help me fix a SyntaxError in my f-strings that use nested quotes
Implement the histogram and main function
The histogram function gives a visual shape to the return distribution right in the terminal. It uses collections.Counter to bucket the data into bins and renders them as ASCII bar charts.
- Add the _print_histogram() function directly below print_report():
def _print_histogram(data: Samples, bins: int = 12, width: int = 35) -> None:
min_val, max_val = min(data), max(data)
bin_width = (max_val - min_val) / bins
bin_counts: collections.Counter = collections.Counter()
for value in data:
bin_idx = min(int((value - min_val) / bin_width), bins - 1)
bin_counts[bin_idx] += 1
max_count = max(bin_counts.values()) if bin_counts else 1
for i in range(bins):
bin_start = min_val + i * bin_width
count = bin_counts[i]
bar_len = int(count / max_count * width)
bar = "\u2588" * bar_len
print(f" {bin_start:>7.1%} | {bar} ({count})")What does this code do?
Each data point is assigned to a bin based on its position between the minimum and maximum values. collections.Counter accumulates the count per bin.
The bar length is normalized so the longest bar fills the available width. The leading underscore in _print_histogram signals it is a private helper, not part of the public API.
Now for the main event. The main() function defines the three investment strategies, runs the optimizer inside a context manager, prints the full report, and demonstrates several Python features you built earlier.
- Add the main execution section below _print_histogram():
# ============================================================
# MAIN EXECUTION
# ============================================================
def main() -> None:
strategies = [
Strategy(
name="Early-Stage",
distribution=LogNormalReturn(mu=0.15, sigma=0.6, label="EarlyStage"),
risk_level=RiskLevel.EXTREME,
),
Strategy(
name="Growth",
distribution=NormalReturn(mu=0.12, sigma=0.25, label="Growth"),
risk_level=RiskLevel.HIGH,
),
Strategy(
name="RealEstate",
distribution=TriangularReturn(low=-0.05, mode=0.08, high=0.20, label="RE"),
risk_level=RiskLevel.MEDIUM,
),
]
print("\n\U0001f4cb Investment Strategies:")
print(f"{'\u2500'*60}")
for s in sorted(strategies):
print(f" {s.name:<15} | {s.distribution.describe()}")
print(f" {'':15} | Risk: {s.risk_level.value} (\u00d7{s.risk_level.multiplier})")
print(f" {'':15} | E[R]: {s.distribution.expected_value():.2%}")
print()
num_trials = 10_000
with SimulationContext("Portfolio Optimization", seed=42, num_trials=num_trials) as sim:
optimal = find_optimal_allocation(
strategies,
num_trials=num_trials,
step=0.1,
objective="sharpe",
seed=sim.seed,
)
sim.results = optimal.returns
print_report(optimal, strategies, fund_size=10_000_000)
print("\n\U0001f3ed Distribution Factory Demo (singledispatch):")
d1 = create_distribution({"type": "normal", "mu": 0.10, "sigma": 0.15})
d2 = create_distribution((-0.05, 0.08, 0.25))
print(f" From dict: {d1}")
print(f" From tuple: {d2}")
print("\n\U0001f9ee Distribution Arithmetic (dunder methods):")
a = NormalReturn(0.10, 0.15, "A")
b = NormalReturn(0.08, 0.20, "B")
combined = a + b
scaled = 2.0 * a
print(f" {a} + {b} = {combined}")
print(f" 2.0 \u00d7 {a} = {scaled}")
print("\n\U0001f4c8 Running Statistics (first 10 trials):")
sample_returns = list(simulate_portfolio(
strategies,
{"Early-Stage": 0.3, "Growth": 0.4, "RealEstate": 0.3},
10,
seed=99,
))
for trial, mean, min_val, max_val in compute_running_stats(iter(sample_returns)):
print(f" Trial {trial:>2}: mean={mean:>7.2%} min={min_val:>7.2%} max={max_val:>7.2%}")
if __name__ == "__main__":
main()What does the main function do?
- Defines three strategies with different distribution types (LogNormal for high-risk early-stage, Normal for growth equity, Triangular for bounded real estate returns).
- Runs the optimizer inside a SimulationContext for reproducibility and clean output formatting.
- Demonstrates create_distribution() singledispatch (creating distributions from dicts and tuples) and operator overloading (adding and scaling NormalReturn objects).
- Shows streaming statistics with compute_running_stats() on a small 10-trial sample.
- Save monte_carlo_engine.py.
Seeing a NameError for a function or class?
Make sure all four parts of the file are in order: Part 1 (distributions), Part 2 (simulation engine), Part 3 (decision analysis), Part 4 (reporting). If main() references a function defined below it, Python will throw a NameError at runtime.
help me fix a NameError in my Monte Carlo engine
✔️ Awesome, I've got everything!
Great. Double check you have saved your file before running it in the next substep.
ⓧ I'd like to double check the full code
Here is the complete monte_carlo_engine.py file. Compare yours against it:
"""
Monte Carlo Decision Analysis Engine
Pure Python Standard Library Implementation
"""
from __future__ import annotations
import random
import statistics
import itertools
import functools
import contextlib
import collections
import enum
import math
import time
from dataclasses import dataclass, field
from typing import Protocol, TypeVar, Iterator, Callable
from abc import ABC, abstractmethod
# ============================================================
# PART 1: Distribution Modeling Layer
# ============================================================
class DistributionError(Exception):
pass
class InvalidParameterError(DistributionError):
def __init__(self, param: str, value: float, reason: str) -> None:
self.param = param
self.value = value
super().__init__(f"Invalid {param}={value}: {reason}")
class RiskLevel(enum.Enum):
LOW = "low"
MEDIUM = "medium"
HIGH = "high"
EXTREME = "extreme"
@property
def multiplier(self) -> float:
multipliers = {
RiskLevel.LOW: 0.5,
RiskLevel.MEDIUM: 1.0,
RiskLevel.HIGH: 1.5,
RiskLevel.EXTREME: 2.0,
}
return multipliers[self]
T = TypeVar("T", bound="Distribution")
class Samplable(Protocol):
def sample(self, n: int, *, seed: int | None = None) -> Iterator[float]: ...
type Samples = list[float]
type Parameters = dict[str, float]
class Distribution(ABC):
@abstractmethod
def sample_one(self, rng: random.Random) -> float: ...
@abstractmethod
def expected_value(self) -> float: ...
@abstractmethod
def describe(self) -> str: ...
def sample(self, n: int, *, seed: int | None = None) -> Iterator[float]:
rng = random.Random(seed)
for _ in itertools.repeat(None, n):
yield self.sample_one(rng)
def __repr__(self) -> str:
return self.describe()
@dataclass(frozen=True, slots=True)
class NormalReturn(Distribution):
mu: float
sigma: float
label: str = field(default="Normal", repr=False)
def __post_init__(self) -> None:
if self.sigma <= 0:
raise InvalidParameterError("sigma", self.sigma, "must be positive")
def sample_one(self, rng: random.Random) -> float:
return rng.gauss(self.mu, self.sigma)
def expected_value(self) -> float:
return self.mu
def describe(self) -> str:
return f"{self.label}(\u03bc={self.mu:.2%}, \u03c3={self.sigma:.2%})"
def __add__(self, other: NormalReturn) -> NormalReturn:
combined_mu = self.mu + other.mu
combined_sigma = math.sqrt(self.sigma**2 + other.sigma**2)
return NormalReturn(combined_mu, combined_sigma, f"{self.label}+{other.label}")
def __mul__(self, scalar: float) -> NormalReturn:
return NormalReturn(self.mu * scalar, self.sigma * abs(scalar), self.label)
def __rmul__(self, scalar: float) -> NormalReturn:
return self.__mul__(scalar)
@dataclass(frozen=True, slots=True)
class TriangularReturn(Distribution):
low: float
mode: float
high: float
label: str = field(default="Triangular", repr=False)
def __post_init__(self) -> None:
if not (self.low <= self.mode <= self.high):
raise InvalidParameterError("mode", self.mode, f"must satisfy low({self.low}) <= mode <= high({self.high})")
def sample_one(self, rng: random.Random) -> float:
return rng.triangular(self.low, self.high, self.mode)
def expected_value(self) -> float:
return (self.low + self.mode + self.high) / 3
def describe(self) -> str:
return f"{self.label}(low={self.low:.2%}, mode={self.mode:.2%}, high={self.high:.2%})"
@dataclass(frozen=True, slots=True)
class LogNormalReturn(Distribution):
mu: float
sigma: float
label: str = field(default="LogNormal", repr=False)
def __post_init__(self) -> None:
if self.sigma <= 0:
raise InvalidParameterError("sigma", self.sigma, "must be positive")
def sample_one(self, rng: random.Random) -> float:
return math.exp(rng.gauss(self.mu, self.sigma)) - 1
def expected_value(self) -> float:
return math.exp(self.mu + self.sigma**2 / 2) - 1
def describe(self) -> str:
return f"{self.label}(\u03bc={self.mu:.3f}, \u03c3={self.sigma:.3f})"
@functools.singledispatch
def create_distribution(spec) -> Distribution:
raise TypeError(f"Cannot create distribution from {type(spec).__name__}")
@create_distribution.register(dict)
def _from_dict(spec: dict) -> Distribution:
match spec:
case {"type": "normal", "mu": mu, "sigma": sigma, **rest}:
return NormalReturn(mu, sigma, rest.get("label", "Normal"))
case {"type": "triangular", "low": low, "mode": mode, "high": high, **rest}:
return TriangularReturn(low, mode, high, rest.get("label", "Triangular"))
case {"type": "lognormal", "mu": mu, "sigma": sigma, **rest}:
return LogNormalReturn(mu, sigma, rest.get("label", "LogNormal"))
case _:
raise DistributionError(f"Unknown distribution spec: {spec}")
@create_distribution.register(tuple)
def _from_tuple(spec: tuple) -> Distribution:
if len(spec) == 3:
return TriangularReturn(*spec)
raise DistributionError(f"Tuple must have 3 elements, got {len(spec)}")
# ============================================================
# PART 2: Simulation Engine
# ============================================================
@dataclass(order=True, kw_only=True)
class Strategy:
sort_index: float = field(init=False, repr=False)
name: str
distribution: Distribution = field(compare=False)
risk_level: RiskLevel = field(compare=False)
allocation: float = field(default=0.0, compare=False)
def __post_init__(self) -> None:
self.sort_index = -self.distribution.expected_value()
@property
def allocation_pct(self) -> str:
return f"{self.allocation:.1%}"
def timer(func: Callable) -> Callable:
@functools.wraps(func)
def wrapper(*args, **kwargs):
start = time.perf_counter()
result = func(*args, **kwargs)
elapsed = time.perf_counter() - start
print(f" \u23f1 {func.__name__} completed in {elapsed:.3f}s")
return result
return wrapper
def validate_allocation(func: Callable) -> Callable:
@functools.wraps(func)
def wrapper(allocations: dict[str, float], *args, **kwargs):
total = sum(allocations.values())
if not math.isclose(total, 1.0, rel_tol=1e-6):
raise ValueError(f"Allocations must sum to 1.0, got {total:.6f}")
return func(allocations, *args, **kwargs)
return wrapper
def retry(max_attempts: int = 3, exceptions: tuple = (Exception,)):
def decorator(func: Callable) -> Callable:
@functools.wraps(func)
def wrapper(*args, **kwargs):
for attempt in range(1, max_attempts + 1):
try:
return func(*args, **kwargs)
except exceptions as e:
if attempt == max_attempts:
raise
print(f" \u26a0 Attempt {attempt} failed: {e}. Retrying...")
return None
return wrapper
return decorator
class SimulationContext:
def __init__(self, name: str, seed: int, num_trials: int) -> None:
self.name = name
self.seed = seed
self.num_trials = num_trials
self.rng: random.Random | None = None
self.results: list[float] = []
def __enter__(self) -> SimulationContext:
self.rng = random.Random(self.seed)
self.results = []
print(f"\n{'='*60}")
print(f" Simulation: {self.name}")
print(f" Trials: {self.num_trials:,} | Seed: {self.seed}")
print(f"{'='*60}")
return self
def __exit__(self, exc_type, exc_val, exc_tb) -> bool:
if exc_type is None:
print(f" \u2713 Completed {len(self.results):,} trials successfully")
else:
print(f" \u2717 Failed: {exc_val}")
return False
@contextlib.contextmanager
def progress_tracker(total: int, label: str = "Progress"):
state = {"completed": 0, "total": total}
def update(n: int = 1) -> None:
state["completed"] += n
if state["completed"] % max(total // 5, 1) == 0:
pct = state["completed"] / state["total"]
print(f" {label}: {pct:.0%}")
try:
yield update
finally:
print(f" {label}: Done ({state['completed']:,}/{state['total']:,})")
def simulate_portfolio(
strategies: list[Strategy],
allocations: dict[str, float],
num_trials: int,
*,
seed: int = 42,
) -> Iterator[float]:
rng = random.Random(seed)
alloc_pairs = [(s, allocations[s.name]) for s in strategies]
for _ in itertools.repeat(None, num_trials):
trial_return = sum(
alloc * s.distribution.sample_one(rng)
for s, alloc in alloc_pairs
)
yield trial_return
def compute_running_stats(returns: Iterator[float]) -> Iterator[tuple]:
t0, t1, t2 = itertools.tee(returns, 3)
running_mean = map(
lambda pair: pair[0] / pair[1],
zip(itertools.accumulate(t0), itertools.count(1))
)
running_max_iter = itertools.accumulate(t1, func=max)
running_min_iter = itertools.accumulate(t2, func=min)
return zip(itertools.count(1), running_mean, running_min_iter, running_max_iter)
# ============================================================
# PART 3: Decision Analysis Framework
# ============================================================
@dataclass
class PortfolioResult:
allocation: dict[str, float]
returns: Samples
_stats_cache: dict = field(default_factory=dict, repr=False)
@functools.cached_property
def mean_return(self) -> float:
return statistics.fmean(self.returns)
@functools.cached_property
def volatility(self) -> float:
return statistics.stdev(self.returns)
@functools.cached_property
def sharpe_ratio(self) -> float:
if (vol := self.volatility) == 0:
return 0.0
return self.mean_return / vol
@functools.cached_property
def var_95(self) -> float:
sorted_returns = sorted(self.returns)
idx = int(len(sorted_returns) * 0.05)
return sorted_returns[idx]
@functools.cached_property
def cvar_95(self) -> float:
sorted_returns = sorted(self.returns)
idx = int(len(sorted_returns) * 0.05)
tail = sorted_returns[:idx] if idx > 0 else sorted_returns[:1]
return statistics.fmean(tail)
@functools.cached_property
def quantiles(self) -> list[float]:
return statistics.quantiles(self.returns, n=4)
def probability_of_loss(self) -> float:
return sum(1 for r in self.returns if r < 0) / len(self.returns)
def __lt__(self, other: PortfolioResult) -> bool:
return self.sharpe_ratio < other.sharpe_ratio
def __eq__(self, other: object) -> bool:
if not isinstance(other, PortfolioResult):
return NotImplemented
return math.isclose(self.sharpe_ratio, other.sharpe_ratio, rel_tol=1e-6)
def generate_allocations(
strategy_names: list[str],
step: float = 0.1,
) -> Iterator[dict[str, float]]:
n = len(strategy_names)
steps = int(1.0 / step) + 1
possible_values = [round(i * step, 2) for i in range(steps)]
for combo in itertools.product(possible_values, repeat=n):
if math.isclose(sum(combo), 1.0, rel_tol=1e-9):
yield dict(zip(strategy_names, combo))
@timer
def find_optimal_allocation(
strategies: list[Strategy],
num_trials: int = 10_000,
step: float = 0.1,
objective: str = "sharpe",
*,
seed: int = 42,
) -> PortfolioResult:
strategy_names = [s.name for s in strategies]
match objective:
case "sharpe":
key_func = lambda r: r.sharpe_ratio
case "return":
key_func = lambda r: r.mean_return
case "risk":
key_func = lambda r: -r.volatility
case "var":
key_func = lambda r: r.var_95
case str() as unknown:
raise ValueError(f"Unknown objective: {unknown}")
best: PortfolioResult | None = None
tested = 0
for allocation in generate_allocations(strategy_names, step):
returns = list(simulate_portfolio(strategies, allocation, num_trials, seed=seed))
result = PortfolioResult(allocation=allocation, returns=returns)
if best is None or (score := key_func(result)) > key_func(best):
best = result
tested += 1
print(f" Tested {tested} improving allocations")
return best
# ============================================================
# PART 4: Reporting & Analysis
# ============================================================
def correlation_matrix(
strategies: list[Strategy],
n_samples: int = 5000,
seed: int = 42,
) -> dict[tuple[str, str], float]:
samples: dict[str, Samples] = {
s.name: list(s.distribution.sample(n_samples, seed=seed))
for s in strategies
}
correlations: dict[tuple[str, str], float] = {}
for (name_a, data_a), (name_b, data_b) in itertools.combinations(
samples.items(), 2
):
corr = statistics.correlation(data_a, data_b)
correlations[(name_a, name_b)] = corr
correlations[(name_b, name_a)] = corr
for name in samples:
correlations[(name, name)] = 1.0
return correlations
def print_report(
result: PortfolioResult,
strategies: list[Strategy],
fund_size: float = 10_000_000,
) -> None:
print(f"\n{'\u2550'*60}")
print(f" \U0001f4ca OPTIMAL PORTFOLIO ANALYSIS")
print(f"{'\u2550'*60}")
print(f"\n \U0001f4b0 Allocation (Fund: ${fund_size:,.0f}):")
for name, weight in result.allocation.items():
dollar_amount = weight * fund_size
print(f" {name:.<30} {weight:>6.1%} (${dollar_amount:>12,.0f})")
print(f"\n \U0001f4c8 Performance Metrics:")
print(f" Expected Return: {result.mean_return:>8.2%}")
print(f" Volatility (\u03c3): {result.volatility:>8.2%}")
print(f" Sharpe Ratio: {result.sharpe_ratio:>8.3f}")
print(f" Value at Risk (95%): {result.var_95:>8.2%}")
print(f" CVaR (95%): {result.cvar_95:>8.2%}")
print(f" P(Loss): {result.probability_of_loss():>8.2%}")
dist = statistics.NormalDist.from_samples(result.returns)
print(f"\n \U0001f4c9 Return Distribution (fitted Normal):")
print(f" Mean: {dist.mean:>8.2%}")
print(f" Stdev: {dist.stdev:>8.2%}")
print(f" P(return > 20%): {1 - dist.cdf(0.20):>8.2%}")
print(f" P(return < -10%): {dist.cdf(-0.10):>8.2%}")
q1, q2, q3 = result.quantiles
print(f"\n \U0001f4ca Quartiles:")
print(f" Q1 (25th): {q1:>8.2%}")
print(f" Q2 (50th): {q2:>8.2%}")
print(f" Q3 (75th): {q3:>8.2%}")
print(f"\n \U0001f4ca Return Distribution Histogram:")
_print_histogram(result.returns)
correlations = correlation_matrix(strategies)
names = [s.name for s in strategies]
print(f"\n \U0001f517 Correlation Matrix:")
header = " " + "".join(f"{n[:10]:>12}" for n in names)
print(header)
for name_a in names:
row = f" {name_a[:10]:<10}"
for name_b in names:
corr = correlations[(name_a, name_b)]
row += f"{corr:>12.3f}"
print(row)
print(f"\n{'\u2550'*60}")
def _print_histogram(data: Samples, bins: int = 12, width: int = 35) -> None:
min_val, max_val = min(data), max(data)
bin_width = (max_val - min_val) / bins
bin_counts: collections.Counter = collections.Counter()
for value in data:
bin_idx = min(int((value - min_val) / bin_width), bins - 1)
bin_counts[bin_idx] += 1
max_count = max(bin_counts.values()) if bin_counts else 1
for i in range(bins):
bin_start = min_val + i * bin_width
count = bin_counts[i]
bar_len = int(count / max_count * width)
bar = "\u2588" * bar_len
print(f" {bin_start:>7.1%} | {bar} ({count})")
# ============================================================
# MAIN EXECUTION
# ============================================================
def main() -> None:
strategies = [
Strategy(
name="Early-Stage",
distribution=LogNormalReturn(mu=0.15, sigma=0.6, label="EarlyStage"),
risk_level=RiskLevel.EXTREME,
),
Strategy(
name="Growth",
distribution=NormalReturn(mu=0.12, sigma=0.25, label="Growth"),
risk_level=RiskLevel.HIGH,
),
Strategy(
name="RealEstate",
distribution=TriangularReturn(low=-0.05, mode=0.08, high=0.20, label="RE"),
risk_level=RiskLevel.MEDIUM,
),
]
print("\n\U0001f4cb Investment Strategies:")
print(f"{'\u2500'*60}")
for s in sorted(strategies):
print(f" {s.name:<15} | {s.distribution.describe()}")
print(f" {'':15} | Risk: {s.risk_level.value} (\u00d7{s.risk_level.multiplier})")
print(f" {'':15} | E[R]: {s.distribution.expected_value():.2%}")
print()
num_trials = 10_000
with SimulationContext("Portfolio Optimization", seed=42, num_trials=num_trials) as sim:
optimal = find_optimal_allocation(
strategies,
num_trials=num_trials,
step=0.1,
objective="sharpe",
seed=sim.seed,
)
sim.results = optimal.returns
print_report(optimal, strategies, fund_size=10_000_000)
print("\n\U0001f3ed Distribution Factory Demo (singledispatch):")
d1 = create_distribution({"type": "normal", "mu": 0.10, "sigma": 0.15})
d2 = create_distribution((-0.05, 0.08, 0.25))
print(f" From dict: {d1}")
print(f" From tuple: {d2}")
print("\n\U0001f9ee Distribution Arithmetic (dunder methods):")
a = NormalReturn(0.10, 0.15, "A")
b = NormalReturn(0.08, 0.20, "B")
combined = a + b
scaled = 2.0 * a
print(f" {a} + {b} = {combined}")
print(f" 2.0 \u00d7 {a} = {scaled}")
print("\n\U0001f4c8 Running Statistics (first 10 trials):")
sample_returns = list(simulate_portfolio(
strategies,
{"Early-Stage": 0.3, "Growth": 0.4, "RealEstate": 0.3},
10,
seed=99,
))
for trial, mean, min_val, max_val in compute_running_stats(iter(sample_returns)):
print(f" Trial {trial:>2}: mean={mean:>7.2%} min={min_val:>7.2%} max={max_val:>7.2%}")
if __name__ == "__main__":
main()Run the engine and interpret the output
Time to see the full engine in action. The optimizer will test all valid allocation combinations and report the one that maximizes the Sharpe ratio.
- Run the engine by executing this command in your terminal:
python3 monte_carlo_engine.pyWhat should I see?
The script takes approximately 5-15 seconds to run (it simulates 10,000 trials across 66 allocation combinations). You should see output in this order:
- An Investment Strategies table showing the three strategies sorted by expected return.
- A Simulation banner with trial count and seed, followed by a timing line.
- The full portfolio report: optimal allocation percentages, dollar amounts, performance metrics (Sharpe ratio, VaR, CVaR), a fitted normal distribution with tail probabilities, quartiles, an ASCII histogram, and a correlation matrix.
- A singledispatch factory demo and a distribution arithmetic demo.
- Running statistics for 10 sample trials showing mean, min, and max converging.
How do I interpret the key metrics?
- Sharpe Ratio measures return per unit of risk. Higher is better. A value above 0.4 indicates reasonable risk-adjusted performance.
- VaR (95%) is the worst return you would expect 95% of the time. A VaR of -15% means only 5% of scenarios are worse than losing 15%.
- CVaR (95%) is the average loss in that worst 5% of cases. It tells you how bad things get when they do go wrong.
- Correlation Matrix values near 0 mean the strategies move independently. Since these distributions use independent random number generators, correlations should be close to zero.
Script takes too long or crashes?
If the script hangs beyond 30 seconds, your machine may be memory-constrained. You can reduce num_trials from 10_000 to 5_000 inside main() for a faster run (results will be noisier but structurally identical).
help me debug my Monte Carlo engine that is crashing or timing out
Secret mission
Sensitivity Analysis & Risk Parity Optimizer
Your engine finds the optimal allocation, but how sensitive is that result to each strategy's performance? And what if you want to equalize risk instead of maximizing returns? Build a sensitivity analysis using linear regression and a risk parity optimizer that balances volatility contributions across strategies.
Clean Up Your Resources
Clean Up Your Resources
Decide whether to keep your resources running, pause them to come back later, or delete them entirely. This project runs entirely locally with Python's standard library, so there are no ongoing costs.
Resources you used:
- Project directory containing monte_carlo_engine.py on your Desktop.
Keep everything running
No action needed. Choose this if you're still actively experimenting with the engine or want to extend it further.
- Your project folder and monte_carlo_engine.py remain on your Desktop, ready to run anytime.
- Python 3.12+ stays installed on your machine for future projects.
Pause - I'll come back to this later
Since this project is a local Python script with no running processes, there is nothing to pause. Your file sits on disk using minimal storage.
- No background services or servers are running from this project.
- You can return to the project anytime by opening the folder in VS Code and running python3 monte_carlo_engine.py again.
Delete - I don't want to use this again
Remove the project directory entirely to free up disk space.
- Delete the project folder by running this command in your terminal:
macOS
cd ~/Desktop
rm -rf monte-carlo-engineWindows
cd ~\Desktop
rmdir /s /q monte-carlo-engine- Confirm the folder is gone by running:
ls ~/DesktopYou should no longer see the project folder in the output.
Nice Work!
Nice Work!
You just built a complete Monte Carlo Decision Analysis Engine using nothing but Python's standard library. From probability distributions to portfolio optimization, every line of code exercised Python's full language power to solve a real analytical problem.
You've learned how to:
- Build a probability distribution modeling layer with ABCs, frozen dataclasses, operator overloading, and singledispatch factories.
- Implement a streaming simulation engine using generators, itertools, decorators, and context managers.
- Create a combinatoric optimization framework with structural pattern matching, cached properties, and comprehensive statistical reporting using NormalDist, correlation matrices, and running statistics.
- Added sensitivity analysis with linear regression and a risk parity optimizer to compare allocation strategies side by side.
Ready to quiz yourself?