A simple tutorial on using Monte Carlo simulation to model stock market returns.
A Simple Yet Complete Tutorial on Estimating Long-Term Investment Returns
Learning Objectives
By the end of this tutorial, you will be able to:
- Understand the fundamentals of Monte Carlo simulation for financial modeling
- Implement a complete investment return simulation using Python
- Interpret probability distributions and risk metrics for investment decisions
- Create meaningful visualizations to communicate financial uncertainty
- Apply these techniques to your own investment analysis
What We'll Build
We'll create a Monte Carlo simulation that estimates the future value of a $100 investment over 10 years, accounting for market volatility and uncertainty. This approach is widely used by financial advisors, portfolio managers, and individual investors to understand potential outcomes and make informed decisions.
Key Concepts Covered
- Monte Carlo Method: Using random sampling to model complex systems
- Investment Returns: How compound growth works with volatility
- Risk Assessment: Understanding percentiles and confidence intervals
- Data Visualization: Creating meaningful charts for financial analysis
1. Setting Up The Environment
First, let's import the essential libraries we'll need for our simulation and visualization:
1import numpy as np
2import pandas as pd
3import matplotlib.pyplot as plt
4
5# Modern, clean styling
6plt.style.use('default')
7plt.rcParams.update({
8 'figure.facecolor': 'white',
9 'axes.facecolor': 'white',
10 'axes.edgecolor': '#CCCCCC',
11 'axes.linewidth': 0.8,
12 'axes.spines.left': True,
13 'axes.spines.bottom': True,
14 'axes.spines.top': False,
15 'axes.spines.right': False,
16 'xtick.bottom': True,
17 'xtick.top': False,
18 'ytick.left': True,
19 'ytick.right': False,
20 'axes.grid': True,
21 'grid.color': '#E5E5E5',
22 'grid.linewidth': 0.5,
23 'font.size': 11,
24 'axes.titlesize': 14,
25 'axes.labelsize': 12,
26 'figure.dpi': 100
27})
1import numpy as np
2import pandas as pd
3import matplotlib.pyplot as plt
4
5# Modern, clean styling
6plt.style.use('default')
7plt.rcParams.update({
8 'figure.facecolor': 'white',
9 'axes.facecolor': 'white',
10 'axes.edgecolor': '#CCCCCC',
11 'axes.linewidth': 0.8,
12 'axes.spines.left': True,
13 'axes.spines.bottom': True,
14 'axes.spines.top': False,
15 'axes.spines.right': False,
16 'xtick.bottom': True,
17 'xtick.top': False,
18 'ytick.left': True,
19 'ytick.right': False,
20 'axes.grid': True,
21 'grid.color': '#E5E5E5',
22 'grid.linewidth': 0.5,
23 'font.size': 11,
24 'axes.titlesize': 14,
25 'axes.labelsize': 12,
26 'figure.dpi': 100
27})
2. Defining The Simulation Parameters
Understanding the Financial Model
Before diving into the code, let's understand what we're modeling:
- Expected Return (): The average annual return expected from the investment (8%)
- Volatility (): How much the returns vary from year to year (15% standard deviation)
- Time Horizon: How long the investment will be held (10 years)
- Simulation Paths: How many different scenarios will be tested (10,000 iterations)
Why These Numbers?
- 8% expected return: Roughly matches historical stock market averages
- 15% volatility: Typical for a diversified stock portfolio
- 10,000 simulations: Provides statistical confidence in the results
The parameters are defined below:
1# Defining parameters for the investment simulation
2INITIAL_INVESTMENT = 100.00
3YEARS = 10
4MU = 0.08 # Expected annual return (8%)
5SIGMA = 0.15 # Annual return volatility (15%)
6N_ITER = 10_000 # Number of Monte Carlo simulation paths
7
8print("Investment Simulation Parameters:")
9print(f" Initial Investment: ${INITIAL_INVESTMENT:,.2f}")
10print(f" Time Horizon: {YEARS} years")
11print(f" Expected Annual Return: {MU:.1%}")
12print(f" Annual Volatility: {SIGMA:.1%}")
13print(f" Number of Simulations: {N_ITER:,}")
1# Defining parameters for the investment simulation
2INITIAL_INVESTMENT = 100.00
3YEARS = 10
4MU = 0.08 # Expected annual return (8%)
5SIGMA = 0.15 # Annual return volatility (15%)
6N_ITER = 10_000 # Number of Monte Carlo simulation paths
7
8print("Investment Simulation Parameters:")
9print(f" Initial Investment: ${INITIAL_INVESTMENT:,.2f}")
10print(f" Time Horizon: {YEARS} years")
11print(f" Expected Annual Return: {MU:.1%}")
12print(f" Annual Volatility: {SIGMA:.1%}")
13print(f" Number of Simulations: {N_ITER:,}")
Investment Simulation Parameters: Initial Investment: $100.00 Time Horizon: 10 years Expected Annual Return: 8.0% Annual Volatility: 15.0% Number of Simulations: 10,000
3. Generating Random Returns
The Heart of Monte Carlo Simulation
Thousands of possible future scenarios will be generated by randomly sampling investment returns from a normal distribution. This is the core of Monte Carlo simulation:
Key Insight: Annual returns are assumed to follow a normal distribution with:
- Mean = 8% (the expected return)
- Standard deviation = 15% (market volatility)
Understanding the Output Structure
- Rows: Each row represents one possible future scenario (simulation path)
- Columns: Each column represents a year in the 10-year horizon
- Values: Each value is a randomly sampled annual return for that year and scenario
1# Setting up reproducible random number generation for consistent results
2rng = np.random.default_rng(seed=42)
3
4# Generating random annual returns for all simulation paths
5# Shape: (N_ITER, YEARS) = (10,000 scenarios × 10 years)
6annual_returns = rng.normal(loc=MU, scale=SIGMA, size=(N_ITER, YEARS))
7
8print(f"Generated Returns Matrix:")
9print(f" Shape: {annual_returns.shape}")
10print(f" Each row = one possible 10-year future")
11print(f" Each column = returns for a specific year across all scenarios")
12print(f"\nSample of first 5 scenarios (first 5 years):")
13print(pd.DataFrame(annual_returns[:5, :5],
14 columns=[f'Year {i+1}' for i in range(5)],
15 index=[f'Scenario {i+1}' for i in range(5)]).round(3))
1# Setting up reproducible random number generation for consistent results
2rng = np.random.default_rng(seed=42)
3
4# Generating random annual returns for all simulation paths
5# Shape: (N_ITER, YEARS) = (10,000 scenarios × 10 years)
6annual_returns = rng.normal(loc=MU, scale=SIGMA, size=(N_ITER, YEARS))
7
8print(f"Generated Returns Matrix:")
9print(f" Shape: {annual_returns.shape}")
10print(f" Each row = one possible 10-year future")
11print(f" Each column = returns for a specific year across all scenarios")
12print(f"\nSample of first 5 scenarios (first 5 years):")
13print(pd.DataFrame(annual_returns[:5, :5],
14 columns=[f'Year {i+1}' for i in range(5)],
15 index=[f'Scenario {i+1}' for i in range(5)]).round(3))
Generated Returns Matrix: Shape: (10000, 10) Each row = one possible 10-year future Each column = returns for a specific year across all scenarios Sample of first 5 scenarios (first 5 years): Year 1 Year 2 Year 3 Year 4 Year 5 Scenario 1 0.126 -0.076 0.193 0.221 -0.213 Scenario 2 0.212 0.197 0.090 0.249 0.150 Scenario 3 0.052 -0.022 0.263 0.057 0.016 Scenario 4 0.401 0.019 0.003 -0.042 0.172 Scenario 5 0.191 0.161 -0.020 0.115 0.098
4. Computing Final Portfolio Values
The Compound Growth Formula
The fundamental principle of compound growth is applied to calculate how the investment grows over time. The mathematical formula is:
Where:
- (capital Pi) means "product of" - multiply all terms together
- is the return in year
- converts a return percentage to a growth factor
Why This Works
- A 10% return means money grows by a factor of 1.10
- A -5% loss means money is multiplied by 0.95
- Over multiple years, all these factors are multiplied together
Example Calculation
If returns of [8%, -2%, 15%] occur over 3 years:
- Growth factors: [1.08, 0.98, 1.15]
- Total growth: 1.08 × 0.98 × 1.15 = 1.217
- 100 × 1.217 = $121.70
1# Converting annual returns to growth factors (1 + return)
2growth_factors_matrix = 1 + annual_returns
3
4# Calculating cumulative growth by multiplying across years (axis=1)
5# This gives the total growth factor for each 10-year scenario
6cumulative_growth_factors = growth_factors_matrix.prod(axis=1)
7
8# Applying compound growth to get final portfolio values
9final_portfolio_values = INITIAL_INVESTMENT * cumulative_growth_factors
10
11print("Example Compound Growth Calculations:")
12print(f" Scenario 1 returns: {annual_returns[0, :].round(3)}")
13print(f" Growth factors: {growth_factors_matrix[0, :].round(3)}")
14print(f" Cumulative growth: {cumulative_growth_factors[0]:.3f}")
15print(f" Final value: ${final_portfolio_values[0]:.2f}")
16print(f"\nSummary Statistics:")
17print(f" Number of scenarios calculated: {len(final_portfolio_values):,}")
18print(f" Minimum final value: ${final_portfolio_values.min():.2f}")
19print(f" Maximum final value: ${final_portfolio_values.max():.2f}")
20print(f" Average final value: ${final_portfolio_values.mean():.2f}")
1# Converting annual returns to growth factors (1 + return)
2growth_factors_matrix = 1 + annual_returns
3
4# Calculating cumulative growth by multiplying across years (axis=1)
5# This gives the total growth factor for each 10-year scenario
6cumulative_growth_factors = growth_factors_matrix.prod(axis=1)
7
8# Applying compound growth to get final portfolio values
9final_portfolio_values = INITIAL_INVESTMENT * cumulative_growth_factors
10
11print("Example Compound Growth Calculations:")
12print(f" Scenario 1 returns: {annual_returns[0, :].round(3)}")
13print(f" Growth factors: {growth_factors_matrix[0, :].round(3)}")
14print(f" Cumulative growth: {cumulative_growth_factors[0]:.3f}")
15print(f" Final value: ${final_portfolio_values[0]:.2f}")
16print(f"\nSummary Statistics:")
17print(f" Number of scenarios calculated: {len(final_portfolio_values):,}")
18print(f" Minimum final value: ${final_portfolio_values.min():.2f}")
19print(f" Maximum final value: ${final_portfolio_values.max():.2f}")
20print(f" Average final value: ${final_portfolio_values.mean():.2f}")
Example Compound Growth Calculations: Scenario 1 returns: [ 0.126 -0.076 0.193 0.221 -0.213 -0.115 0.099 0.033 0.077 -0.048] Growth factors: [1.126 0.924 1.193 1.221 0.787 0.885 1.099 1.033 1.077 0.952] Cumulative growth: 1.228 Final value: $122.84 Summary Statistics: Number of scenarios calculated: 10,000 Minimum final value: $33.48 Maximum final value: $1165.34 Average final value: $214.67
1# Creating comprehensive statistical summary
2portfolio_summary = (pd.Series(final_portfolio_values, name='Final Portfolio Value ($)')
3 .describe(percentiles=[.05, .25, .5, .75, .95])
4 .round(2))
5
6# Reordering for better presentation and adding meaningful labels
7risk_metrics = portfolio_summary.loc[['min', '5%', '25%', '50%', 'mean', '75%', '95%', 'max']]
8risk_metrics.index = ['Worst Case', '5th Percentile (VaR)', '25th Percentile',
9 'Median', 'Expected Value', '75th Percentile',
10 '95th Percentile', 'Best Case']
11
12# Calculating additional useful metrics
13prob_loss = (final_portfolio_values < INITIAL_INVESTMENT).mean() * 100
14prob_double = (final_portfolio_values >= 2 * INITIAL_INVESTMENT).mean() * 100
15expected_return_annualized = (risk_metrics['Expected Value'] / INITIAL_INVESTMENT) ** (1/YEARS) - 1
16
17print("Portfolio Value Distribution After 10 Years")
18print("=" * 50)
19print(risk_metrics.to_string())
20print(f"\nRisk Analysis:")
21print(f" Probability of losing money: {prob_loss:.1f}%")
22print(f" Probability of doubling investment: {prob_double:.1f}%")
23print(f" Expected annualized return: {expected_return_annualized:.1%}")
24print(f" Range of outcomes: ${risk_metrics['Worst Case']:.2f} - ${risk_metrics['Best Case']:.2f}")
1# Creating comprehensive statistical summary
2portfolio_summary = (pd.Series(final_portfolio_values, name='Final Portfolio Value ($)')
3 .describe(percentiles=[.05, .25, .5, .75, .95])
4 .round(2))
5
6# Reordering for better presentation and adding meaningful labels
7risk_metrics = portfolio_summary.loc[['min', '5%', '25%', '50%', 'mean', '75%', '95%', 'max']]
8risk_metrics.index = ['Worst Case', '5th Percentile (VaR)', '25th Percentile',
9 'Median', 'Expected Value', '75th Percentile',
10 '95th Percentile', 'Best Case']
11
12# Calculating additional useful metrics
13prob_loss = (final_portfolio_values < INITIAL_INVESTMENT).mean() * 100
14prob_double = (final_portfolio_values >= 2 * INITIAL_INVESTMENT).mean() * 100
15expected_return_annualized = (risk_metrics['Expected Value'] / INITIAL_INVESTMENT) ** (1/YEARS) - 1
16
17print("Portfolio Value Distribution After 10 Years")
18print("=" * 50)
19print(risk_metrics.to_string())
20print(f"\nRisk Analysis:")
21print(f" Probability of losing money: {prob_loss:.1f}%")
22print(f" Probability of doubling investment: {prob_double:.1f}%")
23print(f" Expected annualized return: {expected_return_annualized:.1%}")
24print(f" Range of outcomes: ${risk_metrics['Worst Case']:.2f} - ${risk_metrics['Best Case']:.2f}")
Portfolio Value Distribution After 10 Years ================================================== Worst Case 33.48 5th Percentile (VaR) 91.25 25th Percentile 144.07 Median 195.04 Expected Value 214.67 75th Percentile 264.35 95th Percentile 404.81 Best Case 1165.34 Risk Analysis: Probability of losing money: 7.6% Probability of doubling investment: 48.0% Expected annualized return: 7.9% Range of outcomes: $33.48 - $1165.34
5. Statistical Analysis and Risk Assessment
Understanding Percentiles and Risk Metrics
The beauty of Monte Carlo simulation lies in its ability to quantify uncertainty. Instead of a single "expected" outcome, a full distribution of possibilities is obtained. The key statistics can be analyzed as follows:
Key Percentiles Explained:
- 5th percentile: Only 5% of scenarios do worse than this (downside risk)
- 25th percentile: First quartile - represents poor but not catastrophic outcomes
- 50th percentile (median): Half of scenarios do better, half do worse
- 75th percentile: Third quartile - represents good outcomes
- 95th percentile: Only 5% of scenarios do better than this (upside potential)
These risk metrics are calculated below:
Understanding Median vs Expected Return
Key Insight: Notice that the median ($195) is lower than the expected value ($215). This is not an error - it's a fundamental characteristic of investment returns.
Why This Happens:
-
Skewed Distribution: Investment returns exhibit positive skewness - there are occasional very large gains that pull the average upward, but losses are bounded (you can't lose more than 100%).
-
Arithmetic vs Geometric: The expected value uses arithmetic averaging of outcomes, while actual compound growth follows geometric progression. A few extremely successful scenarios significantly raise the arithmetic mean.
-
Practical Implication: The median represents the "typical" outcome - half of all scenarios do better, half do worse. The expected value is mathematically correct but influenced by extreme positive outcomes.
Real-World Meaning: If you ran this investment 100 times, you'd be more likely to end up near the median ($195) than the expected value ($215). The expected value includes the impact of those rare scenarios where your portfolio might grow to $500+ or even $1000+.
6. Creating Meaningful Visualizations
We'll create a series of focused visualizations to understand different aspects of our Monte Carlo simulation results. Each chart reveals different insights about the investment risk and return profile.
6.1 Portfolio Value Distribution
The histogram shows the spread of possible outcomes from our Monte Carlo simulation. This visualization helps us understand the likelihood of different portfolio values after 10 years.
1# Portfolio Value Distribution Histogram
2fig, ax = plt.subplots(figsize=(12, 8))
3
4# Filter out values over $600 to focus on the main distribution
5filtered_values = final_portfolio_values[final_portfolio_values <= 600]
6
7# Creating histogram with modern aesthetics
8n, bins, patches = ax.hist(filtered_values, bins=60, alpha=0.7,
9 edgecolor='white', linewidth=0.5, color='#3498db')
10
11# Adding vertical lines for key percentiles with better positioning
12percentiles = [5, 25, 50, 75, 95]
13percentile_values = np.percentile(final_portfolio_values, percentiles)
14colors = ['#e74c3c', '#f39c12', '#27ae60', '#f39c12', '#e74c3c']
15labels = ['5th%', '25th%', 'Median', '75th%', '95th%']
16
17max_height = n.max()
18for i, (perc, val, color, label) in enumerate(zip(percentiles, percentile_values, colors, labels)):
19 if val <= 600:
20 ax.axvline(val, color=color, linestyle='-', linewidth=2, alpha=0.8)
21 ax.text(val, max_height * (0.85 - i * 0.1), f'{label}\n${val:.0f}',
22 ha='center', va='top', fontweight='600', fontsize=10,
23 bbox=dict(boxstyle="round,pad=0.4", facecolor=color, alpha=0.15, edgecolor=color))
24
25# Adding initial investment line for reference
26ax.axvline(INITIAL_INVESTMENT, color='#2c3e50', linestyle='-', linewidth=2)
27ax.text(INITIAL_INVESTMENT, max_height * 0.95, f'Initial\n${INITIAL_INVESTMENT:.0f}',
28 ha='center', va='top', fontweight='600', fontsize=10,
29 bbox=dict(boxstyle="round,pad=0.4", facecolor='#95a5a6', alpha=0.15, edgecolor='#2c3e50'))
30
31# Calculate excluded scenarios
32excluded_count = len(final_portfolio_values) - len(filtered_values)
33excluded_pct = (excluded_count / len(final_portfolio_values)) * 100
34
35ax.set_title(f"Portfolio Value Distribution After {YEARS} Years\n"
36 f"Expected Return: {MU:.1%} • Volatility: {SIGMA:.1%} • {N_ITER:,} Simulations",
37 fontsize=16, fontweight='600', pad=20, color='#2c3e50')
38ax.set_xlabel("Final Portfolio Value ($)", fontsize=13, fontweight='500', color='#2c3e50')
39ax.set_ylabel("Number of Simulations", fontsize=13, fontweight='500', color='#2c3e50')
40
41# Adding summary statistics in a cleaner box
42stats_text = f"Expected Value: ${final_portfolio_values.mean():.0f}\n" \
43 f"Median: ${np.median(final_portfolio_values):.0f}\n" \
44 f"Loss Probability: {prob_loss:.1f}%\n" \
45 f"Double Probability: {prob_double:.1f}%"
46
47ax.text(0.97, 0.97, stats_text, transform=ax.transAxes, fontsize=11, fontweight='500',
48 va='top', ha='right', color='#2c3e50',
49 bbox=dict(boxstyle="round,pad=0.5", facecolor='white', alpha=0.9,
50 edgecolor='#bdc3c7', linewidth=1))
51
52# Clean up the plot
53ax.grid(True, alpha=0.3, color='#ecf0f1')
54ax.set_axisbelow(True)
55
56plt.tight_layout()
57plt.show()
1# Portfolio Value Distribution Histogram
2fig, ax = plt.subplots(figsize=(12, 8))
3
4# Filter out values over $600 to focus on the main distribution
5filtered_values = final_portfolio_values[final_portfolio_values <= 600]
6
7# Creating histogram with modern aesthetics
8n, bins, patches = ax.hist(filtered_values, bins=60, alpha=0.7,
9 edgecolor='white', linewidth=0.5, color='#3498db')
10
11# Adding vertical lines for key percentiles with better positioning
12percentiles = [5, 25, 50, 75, 95]
13percentile_values = np.percentile(final_portfolio_values, percentiles)
14colors = ['#e74c3c', '#f39c12', '#27ae60', '#f39c12', '#e74c3c']
15labels = ['5th%', '25th%', 'Median', '75th%', '95th%']
16
17max_height = n.max()
18for i, (perc, val, color, label) in enumerate(zip(percentiles, percentile_values, colors, labels)):
19 if val <= 600:
20 ax.axvline(val, color=color, linestyle='-', linewidth=2, alpha=0.8)
21 ax.text(val, max_height * (0.85 - i * 0.1), f'{label}\n${val:.0f}',
22 ha='center', va='top', fontweight='600', fontsize=10,
23 bbox=dict(boxstyle="round,pad=0.4", facecolor=color, alpha=0.15, edgecolor=color))
24
25# Adding initial investment line for reference
26ax.axvline(INITIAL_INVESTMENT, color='#2c3e50', linestyle='-', linewidth=2)
27ax.text(INITIAL_INVESTMENT, max_height * 0.95, f'Initial\n${INITIAL_INVESTMENT:.0f}',
28 ha='center', va='top', fontweight='600', fontsize=10,
29 bbox=dict(boxstyle="round,pad=0.4", facecolor='#95a5a6', alpha=0.15, edgecolor='#2c3e50'))
30
31# Calculate excluded scenarios
32excluded_count = len(final_portfolio_values) - len(filtered_values)
33excluded_pct = (excluded_count / len(final_portfolio_values)) * 100
34
35ax.set_title(f"Portfolio Value Distribution After {YEARS} Years\n"
36 f"Expected Return: {MU:.1%} • Volatility: {SIGMA:.1%} • {N_ITER:,} Simulations",
37 fontsize=16, fontweight='600', pad=20, color='#2c3e50')
38ax.set_xlabel("Final Portfolio Value ($)", fontsize=13, fontweight='500', color='#2c3e50')
39ax.set_ylabel("Number of Simulations", fontsize=13, fontweight='500', color='#2c3e50')
40
41# Adding summary statistics in a cleaner box
42stats_text = f"Expected Value: ${final_portfolio_values.mean():.0f}\n" \
43 f"Median: ${np.median(final_portfolio_values):.0f}\n" \
44 f"Loss Probability: {prob_loss:.1f}%\n" \
45 f"Double Probability: {prob_double:.1f}%"
46
47ax.text(0.97, 0.97, stats_text, transform=ax.transAxes, fontsize=11, fontweight='500',
48 va='top', ha='right', color='#2c3e50',
49 bbox=dict(boxstyle="round,pad=0.5", facecolor='white', alpha=0.9,
50 edgecolor='#bdc3c7', linewidth=1))
51
52# Clean up the plot
53ax.grid(True, alpha=0.3, color='#ecf0f1')
54ax.set_axisbelow(True)
55
56plt.tight_layout()
57plt.show()
6.2 Box Plot Summary
The box plot provides a concise statistical summary, highlighting the quartiles, median, and outliers in our portfolio value distribution.
1# Box Plot showing percentile ranges with modern styling
2fig, ax = plt.subplots(figsize=(10, 8))
3
4box_data = [final_portfolio_values]
5bp = ax.boxplot(box_data, patch_artist=True, tick_labels=['Portfolio Values'],
6 medianprops=dict(color='#27ae60', linewidth=2.5),
7 boxprops=dict(facecolor='#3498db', alpha=0.6, edgecolor='#2980b9', linewidth=1.5),
8 whiskerprops=dict(color='#34495e', linewidth=1.5),
9 capprops=dict(color='#34495e', linewidth=1.5),
10 flierprops=dict(marker='o', markerfacecolor='#e74c3c', alpha=0.5, markersize=4))
11
12# Adding percentile annotations with cleaner styling
13percentiles_to_show = [5, 25, 50, 75, 95]
14perc_values = np.percentile(final_portfolio_values, percentiles_to_show)
15colors = ['#e74c3c', '#f39c12', '#27ae60', '#f39c12', '#e74c3c']
16
17for i, (perc, val, color) in enumerate(zip(percentiles_to_show, perc_values, colors)):
18 ax.annotate(f'{perc}th: ${val:.0f}',
19 xy=(1, val), xytext=(1.15, val),
20 fontsize=11, ha='left', va='center', fontweight='600',
21 color=color,
22 arrowprops=dict(arrowstyle='->', color=color, alpha=0.7, lw=1.2))
23
24ax.set_title('Portfolio Value Distribution Summary',
25 fontweight='600', fontsize=16, color='#2c3e50', pad=20)
26ax.set_ylabel('Portfolio Value ($)', fontweight='500', fontsize=13, color='#2c3e50')
27ax.grid(True, alpha=0.3, color='#ecf0f1', axis='y')
28ax.set_axisbelow(True)
29
30# Remove top and right spines for cleaner look
31ax.spines['top'].set_visible(False)
32ax.spines['right'].set_visible(False)
33ax.spines['left'].set_color('#bdc3c7')
34ax.spines['bottom'].set_color('#bdc3c7')
35
36plt.tight_layout()
37plt.show()
1# Box Plot showing percentile ranges with modern styling
2fig, ax = plt.subplots(figsize=(10, 8))
3
4box_data = [final_portfolio_values]
5bp = ax.boxplot(box_data, patch_artist=True, tick_labels=['Portfolio Values'],
6 medianprops=dict(color='#27ae60', linewidth=2.5),
7 boxprops=dict(facecolor='#3498db', alpha=0.6, edgecolor='#2980b9', linewidth=1.5),
8 whiskerprops=dict(color='#34495e', linewidth=1.5),
9 capprops=dict(color='#34495e', linewidth=1.5),
10 flierprops=dict(marker='o', markerfacecolor='#e74c3c', alpha=0.5, markersize=4))
11
12# Adding percentile annotations with cleaner styling
13percentiles_to_show = [5, 25, 50, 75, 95]
14perc_values = np.percentile(final_portfolio_values, percentiles_to_show)
15colors = ['#e74c3c', '#f39c12', '#27ae60', '#f39c12', '#e74c3c']
16
17for i, (perc, val, color) in enumerate(zip(percentiles_to_show, perc_values, colors)):
18 ax.annotate(f'{perc}th: ${val:.0f}',
19 xy=(1, val), xytext=(1.15, val),
20 fontsize=11, ha='left', va='center', fontweight='600',
21 color=color,
22 arrowprops=dict(arrowstyle='->', color=color, alpha=0.7, lw=1.2))
23
24ax.set_title('Portfolio Value Distribution Summary',
25 fontweight='600', fontsize=16, color='#2c3e50', pad=20)
26ax.set_ylabel('Portfolio Value ($)', fontweight='500', fontsize=13, color='#2c3e50')
27ax.grid(True, alpha=0.3, color='#ecf0f1', axis='y')
28ax.set_axisbelow(True)
29
30# Remove top and right spines for cleaner look
31ax.spines['top'].set_visible(False)
32ax.spines['right'].set_visible(False)
33ax.spines['left'].set_color('#bdc3c7')
34ax.spines['bottom'].set_color('#bdc3c7')
35
36plt.tight_layout()
37plt.show()
6.3 Cumulative Distribution Function (CDF)
The CDF shows the probability of achieving different portfolio values. This helps answer questions like "What's the probability my portfolio will be worth at least $200?"
1# Cumulative Distribution Function (CDF) with modern styling
2fig, ax = plt.subplots(figsize=(12, 8))
3
4sorted_values = np.sort(final_portfolio_values)
5cumulative_prob = np.arange(1, len(sorted_values) + 1) / len(sorted_values)
6
7# Main CDF line with gradient-like effect
8ax.plot(sorted_values, cumulative_prob * 100, linewidth=3, color='#3498db',
9 label='Cumulative Probability', alpha=0.9)
10
11# Reference lines
12ax.axvline(INITIAL_INVESTMENT, color='#e74c3c', linestyle='--', alpha=0.8,
13 linewidth=2, label='Break-even ($100)')
14ax.axhline(50, color='#95a5a6', linestyle=':', alpha=0.7,
15 linewidth=2, label='50% Probability')
16
17# Key probability markers with modern styling
18key_values = [150, 200, 250, 300]
19marker_colors = ['#e67e22', '#f39c12', '#27ae60', '#8e44ad']
20
21for val, color in zip(key_values, marker_colors):
22 prob = (final_portfolio_values <= val).mean() * 100
23 ax.plot(val, prob, 'o', markersize=8, color=color,
24 markeredgecolor='white', markeredgewidth=2, zorder=5)
25 ax.annotate(f'${val}: {prob:.0f}%',
26 xy=(val, prob), xytext=(val+25, prob+8),
27 fontsize=11, ha='left', fontweight='600', color=color,
28 bbox=dict(boxstyle="round,pad=0.3", facecolor=color, alpha=0.15,
29 edgecolor=color, linewidth=1),
30 arrowprops=dict(arrowstyle='->', color=color, alpha=0.7, lw=1.2))
31
32ax.set_title('Cumulative Probability Distribution',
33 fontweight='600', fontsize=16, color='#2c3e50', pad=20)
34ax.set_xlabel('Portfolio Value ($)', fontweight='500', fontsize=13, color='#2c3e50')
35ax.set_ylabel('Probability (%)', fontweight='500', fontsize=13, color='#2c3e50')
36
37# Modern legend styling
38ax.legend(fontsize=11, frameon=True, fancybox=True, shadow=False,
39 facecolor='white', edgecolor='#bdc3c7', framealpha=0.9)
40
41# Clean grid and spines
42ax.grid(True, alpha=0.3, color='#ecf0f1')
43ax.set_axisbelow(True)
44ax.spines['top'].set_visible(False)
45ax.spines['right'].set_visible(False)
46ax.spines['left'].set_color('#bdc3c7')
47ax.spines['bottom'].set_color('#bdc3c7')
48
49plt.tight_layout()
50plt.show()
1# Cumulative Distribution Function (CDF) with modern styling
2fig, ax = plt.subplots(figsize=(12, 8))
3
4sorted_values = np.sort(final_portfolio_values)
5cumulative_prob = np.arange(1, len(sorted_values) + 1) / len(sorted_values)
6
7# Main CDF line with gradient-like effect
8ax.plot(sorted_values, cumulative_prob * 100, linewidth=3, color='#3498db',
9 label='Cumulative Probability', alpha=0.9)
10
11# Reference lines
12ax.axvline(INITIAL_INVESTMENT, color='#e74c3c', linestyle='--', alpha=0.8,
13 linewidth=2, label='Break-even ($100)')
14ax.axhline(50, color='#95a5a6', linestyle=':', alpha=0.7,
15 linewidth=2, label='50% Probability')
16
17# Key probability markers with modern styling
18key_values = [150, 200, 250, 300]
19marker_colors = ['#e67e22', '#f39c12', '#27ae60', '#8e44ad']
20
21for val, color in zip(key_values, marker_colors):
22 prob = (final_portfolio_values <= val).mean() * 100
23 ax.plot(val, prob, 'o', markersize=8, color=color,
24 markeredgecolor='white', markeredgewidth=2, zorder=5)
25 ax.annotate(f'${val}: {prob:.0f}%',
26 xy=(val, prob), xytext=(val+25, prob+8),
27 fontsize=11, ha='left', fontweight='600', color=color,
28 bbox=dict(boxstyle="round,pad=0.3", facecolor=color, alpha=0.15,
29 edgecolor=color, linewidth=1),
30 arrowprops=dict(arrowstyle='->', color=color, alpha=0.7, lw=1.2))
31
32ax.set_title('Cumulative Probability Distribution',
33 fontweight='600', fontsize=16, color='#2c3e50', pad=20)
34ax.set_xlabel('Portfolio Value ($)', fontweight='500', fontsize=13, color='#2c3e50')
35ax.set_ylabel('Probability (%)', fontweight='500', fontsize=13, color='#2c3e50')
36
37# Modern legend styling
38ax.legend(fontsize=11, frameon=True, fancybox=True, shadow=False,
39 facecolor='white', edgecolor='#bdc3c7', framealpha=0.9)
40
41# Clean grid and spines
42ax.grid(True, alpha=0.3, color='#ecf0f1')
43ax.set_axisbelow(True)
44ax.spines['top'].set_visible(False)
45ax.spines['right'].set_visible(False)
46ax.spines['left'].set_color('#bdc3c7')
47ax.spines['bottom'].set_color('#bdc3c7')
48
49plt.tight_layout()
50plt.show()
6.4 Portfolio Growth Over Time
This visualization shows how portfolio values evolve year by year, demonstrating the compound growth effect and the confidence band around the expected path.
1# Portfolio growth over time with modern styling
2fig, ax = plt.subplots(figsize=(12, 8))
3
4# Calculate portfolio values for each year for sample scenarios
5sample_scenarios = annual_returns[:100, :] # First 100 scenarios
6portfolio_evolution = np.zeros((100, YEARS + 1)) # +1 for initial value
7portfolio_evolution[:, 0] = INITIAL_INVESTMENT
8
9for year in range(YEARS):
10 portfolio_evolution[:, year + 1] = portfolio_evolution[:, year] * (1 + sample_scenarios[:, year])
11
12# Plot sample paths (very light gray lines)
13for i in range(20):
14 ax.plot(np.arange(YEARS + 1), portfolio_evolution[i, :],
15 alpha=0.15, color='#95a5a6', linewidth=0.8, zorder=1)
16
17# Calculate evolution for all scenarios for confidence bands
18all_portfolio_evolution = np.zeros((N_ITER, YEARS + 1))
19all_portfolio_evolution[:, 0] = INITIAL_INVESTMENT
20
21for year in range(YEARS):
22 all_portfolio_evolution[:, year + 1] = all_portfolio_evolution[:, year] * (1 + annual_returns[:, year])
23
24# Add confidence bands with gradient-like colors
25percentile_5 = np.percentile(all_portfolio_evolution, 5, axis=0)
26percentile_95 = np.percentile(all_portfolio_evolution, 95, axis=0)
27percentile_25 = np.percentile(all_portfolio_evolution, 25, axis=0)
28percentile_75 = np.percentile(all_portfolio_evolution, 75, axis=0)
29
30# Fill confidence bands with different opacities
31ax.fill_between(np.arange(YEARS + 1), percentile_5, percentile_95,
32 alpha=0.15, color='#3498db', label='90% Confidence Band', zorder=2)
33ax.fill_between(np.arange(YEARS + 1), percentile_25, percentile_75,
34 alpha=0.25, color='#3498db', label='50% Confidence Band', zorder=3)
35
36# Plot median line
37median_evolution = np.percentile(all_portfolio_evolution, 50, axis=0)
38ax.plot(np.arange(YEARS + 1), median_evolution, linewidth=3, color='#27ae60',
39 linestyle='--', label='Median Growth Path', zorder=4)
40
41# Plot average evolution
42mean_evolution = all_portfolio_evolution.mean(axis=0)
43ax.plot(np.arange(YEARS + 1), mean_evolution, linewidth=3, color='#e74c3c',
44 label='Expected Growth Path', zorder=5)
45
46# Add initial investment reference line
47ax.axhline(INITIAL_INVESTMENT, color='#34495e', linestyle=':', alpha=0.7,
48 linewidth=1.5, label='Initial Investment')
49
50ax.set_title('Portfolio Growth Trajectories Over Time',
51 fontweight='600', fontsize=16, color='#2c3e50', pad=20)
52ax.set_xlabel('Year', fontweight='500', fontsize=13, color='#2c3e50')
53ax.set_ylabel('Portfolio Value ($)', fontweight='500', fontsize=13, color='#2c3e50')
54
55# Modern legend styling
56ax.legend(fontsize=11, frameon=True, fancybox=True, shadow=False,
57 facecolor='white', edgecolor='#bdc3c7', framealpha=0.9,
58 loc='upper left')
59
60# Clean grid and spines
61ax.grid(True, alpha=0.3, color='#ecf0f1')
62ax.set_axisbelow(True)
63ax.spines['top'].set_visible(False)
64ax.spines['right'].set_visible(False)
65ax.spines['left'].set_color('#bdc3c7')
66ax.spines['bottom'].set_color('#bdc3c7')
67
68# Set reasonable y-axis limits to avoid extreme outliers
69ax.set_ylim(0, 600)
70
71plt.tight_layout()
72plt.show()
1# Portfolio growth over time with modern styling
2fig, ax = plt.subplots(figsize=(12, 8))
3
4# Calculate portfolio values for each year for sample scenarios
5sample_scenarios = annual_returns[:100, :] # First 100 scenarios
6portfolio_evolution = np.zeros((100, YEARS + 1)) # +1 for initial value
7portfolio_evolution[:, 0] = INITIAL_INVESTMENT
8
9for year in range(YEARS):
10 portfolio_evolution[:, year + 1] = portfolio_evolution[:, year] * (1 + sample_scenarios[:, year])
11
12# Plot sample paths (very light gray lines)
13for i in range(20):
14 ax.plot(np.arange(YEARS + 1), portfolio_evolution[i, :],
15 alpha=0.15, color='#95a5a6', linewidth=0.8, zorder=1)
16
17# Calculate evolution for all scenarios for confidence bands
18all_portfolio_evolution = np.zeros((N_ITER, YEARS + 1))
19all_portfolio_evolution[:, 0] = INITIAL_INVESTMENT
20
21for year in range(YEARS):
22 all_portfolio_evolution[:, year + 1] = all_portfolio_evolution[:, year] * (1 + annual_returns[:, year])
23
24# Add confidence bands with gradient-like colors
25percentile_5 = np.percentile(all_portfolio_evolution, 5, axis=0)
26percentile_95 = np.percentile(all_portfolio_evolution, 95, axis=0)
27percentile_25 = np.percentile(all_portfolio_evolution, 25, axis=0)
28percentile_75 = np.percentile(all_portfolio_evolution, 75, axis=0)
29
30# Fill confidence bands with different opacities
31ax.fill_between(np.arange(YEARS + 1), percentile_5, percentile_95,
32 alpha=0.15, color='#3498db', label='90% Confidence Band', zorder=2)
33ax.fill_between(np.arange(YEARS + 1), percentile_25, percentile_75,
34 alpha=0.25, color='#3498db', label='50% Confidence Band', zorder=3)
35
36# Plot median line
37median_evolution = np.percentile(all_portfolio_evolution, 50, axis=0)
38ax.plot(np.arange(YEARS + 1), median_evolution, linewidth=3, color='#27ae60',
39 linestyle='--', label='Median Growth Path', zorder=4)
40
41# Plot average evolution
42mean_evolution = all_portfolio_evolution.mean(axis=0)
43ax.plot(np.arange(YEARS + 1), mean_evolution, linewidth=3, color='#e74c3c',
44 label='Expected Growth Path', zorder=5)
45
46# Add initial investment reference line
47ax.axhline(INITIAL_INVESTMENT, color='#34495e', linestyle=':', alpha=0.7,
48 linewidth=1.5, label='Initial Investment')
49
50ax.set_title('Portfolio Growth Trajectories Over Time',
51 fontweight='600', fontsize=16, color='#2c3e50', pad=20)
52ax.set_xlabel('Year', fontweight='500', fontsize=13, color='#2c3e50')
53ax.set_ylabel('Portfolio Value ($)', fontweight='500', fontsize=13, color='#2c3e50')
54
55# Modern legend styling
56ax.legend(fontsize=11, frameon=True, fancybox=True, shadow=False,
57 facecolor='white', edgecolor='#bdc3c7', framealpha=0.9,
58 loc='upper left')
59
60# Clean grid and spines
61ax.grid(True, alpha=0.3, color='#ecf0f1')
62ax.set_axisbelow(True)
63ax.spines['top'].set_visible(False)
64ax.spines['right'].set_visible(False)
65ax.spines['left'].set_color('#bdc3c7')
66ax.spines['bottom'].set_color('#bdc3c7')
67
68# Set reasonable y-axis limits to avoid extreme outliers
69ax.set_ylim(0, 600)
70
71plt.tight_layout()
72plt.show()
7. Conclusion: The Power of Monte Carlo Simulation
Monte Carlo simulation transforms investment uncertainty from guesswork into quantified risk assessment. By running 10,000 possible scenarios, we've mapped the full landscape of potential outcomes for our investment.
Key Monte Carlo Insights:
- Probabilistic Thinking: Rather than a single "expected" return, we now understand the full distribution of possibilities
- Risk Quantification: We can precisely state there's a 7.6% chance of losing money and a 48% chance of doubling our investment
- Confidence Intervals: We're 90% confident our final portfolio will be between 405
Why Monte Carlo Works:
- Captures Uncertainty: Markets are inherently random - Monte Carlo embraces this reality rather than ignoring it
- Compound Effects: Shows how volatility compounds over time, revealing both upside potential and downside risk
- Decision Support: Provides the statistical foundation for rational investment decisions
The Monte Carlo Advantage:
Traditional financial planning might say "expect 8% returns." Monte Carlo simulation reveals that while 8% is the average, actual outcomes range dramatically. This knowledge is power - it enables better risk management, more realistic expectations, and informed decision-making.
Monte Carlo simulation is not just a mathematical exercise; it's a lens for understanding uncertainty in any complex system where randomness plays a crucial role.

About the author: Michael Brenndoerfer
All opinions expressed here are my own and do not reflect the views of my employer.
Michael currently works as an Associate Director of Data Science at EQT Partners in Singapore, where he drives AI and data initiatives across private capital investments.
With over a decade of experience spanning private equity, management consulting, and software engineering, he specializes in building and scaling analytics capabilities from the ground up. He has published research in leading AI conferences and holds expertise in machine learning, natural language processing, and value creation through data.
Related Content

What are AI Agents, Really?
A comprehensive guide to understanding AI agents, their building blocks, and how they differ from agentic workflows and agent swarms.

Understanding the Model Context Protocol (MCP)
A deep dive into how MCP makes tool use with LLMs easier, cleaner, and more standardized.

Why Temperature=0 Doesn't Guarantee Determinism in LLMs
An exploration of why setting temperature to zero doesn't eliminate all randomness in large language model outputs.
Stay Updated
Get notified when I publish new articles on data and AI, private equity, technology, and more.