### Simulating future stock prices using Monte Carlo methods in Python

Using Monte Carlo methods, we'll write a quick simulation to predict future stock price outcomes for Apple ($AAPL) using Python. You can read more about Monte Carlo simulation (in a finance context) here. 1) Pull the data First, we can import the libraries, and pull the historical stock data for Apple. For this example, I picked the last ~10 years, although it would be valuable to test sensitivities of different ranges as this alone is subjective. # Import required libraries import math import matplotlib.pyplot as plt import numpy as np from pandas_datareader import data apple = data.DataReader('AAPL', 'yahoo',start='1/1/2009') apple.head() High Low Open Close Volume Adj Close Date 2009-01-02 13.005714 12.165714 12.268572 12.964286 186503800.0 11.357091 2009-01-05 13.740000 13.244286 13.310000 13.511429 295402100.0 11.836406 2009-01-06 13.881429 13.198571 13.707143 13.288571 322327600.0 11.641174 2009-01-07 13.214286 12.894286 13.115714 13.001429 188262200.0 11.389630 2009-01-08 13.307143 12.862857 12.918571 13.242857 168375200.0 11.601127 #Next, we calculate the number of days that have elapsed in our chosen time window time_elapsed = (apple.index[-1] - apple.index[0]).days 2) Calculate the compounded annualized growth rate over the length of the dataset + standard deviation (to feed into simulation) #Current price / first record (e.g. price at beginning of 2009) #provides us with the total growth % total_growth = (apple['Adj Close'][-1] / apple['Adj Close'][1]) #Next, we want to annualize this percentage #First, we convert our time elapsed to the # of years elapsed number_of_years = time_elapsed / 365.0 #Second, we can raise the total growth to the inverse of the # of years #(e.g. ~1/10 at time of writing) to annualize our growth rate cagr = total_growth ** (1/number_of_years) - 1 #Now that we have the mean annual growth rate above, #we'll also need to calculate the standard deviation of the #daily price changes std_dev = apple['Adj Close'].pct_change().std() #Next, because there are roughy ~252 trading days in a year, #we'll need to scale this by an annualization factor #reference: https://www.fool.com/knowledge-center/how-to-calculate-annualized-volatility.aspx number_of_trading_days = 252 std_dev = std_dev * math.sqrt(number_of_trading_days) #From here, we have our two inputs needed to generate random #values in our simulation print ("cagr (mean returns) : ", str(round(cagr,4))) print ("std_dev (standard deviation of return : )", str(round(std_dev,4))) cagr (mean returns) : 0.3113 std_dev (standard deviation of return : ) 0.2671 3) Generate random values, run Monte Carlo simulation #Generate random values for 1 year's worth of trading (252 days), #using numpy and assuming a normal distribution daily_return_percentages = np.random.normal(cagr/number_of_trading_days, std_dev/math.sqrt(number_of_trading_days) number_of_trading_days)+1 #Now that we have created a random series of future #daily return %s, we can simply apply these forward-looking #to our last stock price in the window, effectively carrying forward #a price prediction for the next year #This distribution is known as a 'random walk' price_series = [apple['Adj Close'][-1]] for j in daily_return_percentages: price_series.append(price_series[-1] * j) #Great, now we can plot of single 'random walk' of stock prices plt.plot(price_series) plt.show() #Now that we've created a single random walk above, #we can simulate this process over a large sample size to #get a better sense of the true expected distribution number_of_trials = 3000 #set up an additional array to collect all possible #closing prices in last day of window. #We can toss this into a histogram #to get a clearer sense of possible outcomes closing_prices = [] for i in range(number_of_trials): #calculate randomized return percentages following our normal distribution #and using the mean / std dev we calculated above daily_return_percentages = np.random.normal(cagr/number_of_trading_days, std_dev/math.sqrt(number_of_trading_days), number_of_trading_days)+1 price_series = [apple['Adj Close'][-1]] for j in daily_return_percentages: #extrapolate price out for next year price_series.append(price_series[-1] * j) #append closing prices in last day of window for histogram closing_prices.append(price_series[-1]) #plot all random walks plt.plot(price_series) plt.show() #plot histogram plt.hist(closing_prices,bins=40) plt.show() 4) Analyze results #from here, we can check the mean of all ending prices #allowing us to arrive at the most probable ending point mean_end_price = round(np.mean(closing_prices),2) print("Expected price: ", str(mean_end_price)) Expected price: 276.05 #lastly, we can split the distribution into percentiles #to help us gauge risk vs. reward #Pull top 10% of possible outcomes top_ten = np.percentile(closing_prices,100-10) #Pull bottom 10% of possible outcomes bottom_ten = np.percentile(closing_prices,10); #create histogram again plt.hist(closing_prices,bins=40) #append w/ top 10% line plt.axvline(top_ten,color='r',linestyle='dashed',linewidth=2) #append w/ bottom 10% line plt.axvline(bottom_ten,color='r',linestyle='dashed',linewidth=2) #append with current price plt.axhline(apple['Adj Close'][-1],color='g', linestyle='dashed',linewidth=2) plt.show() Here, we can see (based solely on using Monte Carlo simulation, of course) there looks to be more upside than downside for the next year, with the expected price running about$270 and only a 10 percent chance of the price landing below 189.

There are of course many, many things that drive a stock price beyond historical percent changes (meaning I would not make investment decisions solely off Monte Carlo simulation!), but this offers a thorough example of using Monte Carlo simulation to better understand the distribution of possible outcomes.

### Ace your next data science interview

Get better at data science interviews by solving a few questions per week