<jupyter_start><jupyter_text>Stock prices simulation:<jupyter_code>def stock_price_generator (T, n ,m, r, S0, k, V0, sigma, theta, kappa, rho, generator):
    dt = T / n
    
    # Brownian motions:
    dw_v = generator.normal(size=(m, n)) * np.sqrt(dt)
    dw_i = generator.normal(size=(m, n)) * np.sqrt(dt)
    #dw_v = np.random.normal(size=(m, n)) * np.sqrt(dt)
    #dw_i = np.random.normal(size=(m, n)) * np.sqrt(dt)

    dw_s = rho * dw_v + np.sqrt(1.0 - rho ** 2) * dw_i

    # Perform time evolution 
    s = np.empty((m, n + 1)) # initialisation stock prices vector
    s[:, 0] = S0

    v = np.ones(m) * V0

    for t in range(n):
        dv = kappa * (theta - v) * dt + sigma * np.sqrt(v) * dw_v[:, t]
        ds = r * s[:, t] * dt + np.sqrt(v) * s[:, t] * dw_s[:, t]

        v = np.clip(v + dv, a_min=0.0, a_max=None)
        s[:, t + 1] = s[:, t] + ds
        
        
    return s[:,-1]
    
#s = stock_price_generator (, n ,m, r, S0, k, V0, sigma, theta, kappa, rho)<jupyter_output><empty_output><jupyter_text>Plots stock prices simulations:<jupyter_code>#plot_mc (s,,n)<jupyter_output><empty_output><jupyter_text>Payoff:<jupyter_code># function which finds the expected payoff 
def find_expected_payoff(stock_paths, k, r, T):
    final = stock_paths - k
    payoff = [max(x,0) for x in final] # one payoff for each simulation
    payoff = np.asarray(payoff)
    expected_payoff = payoff.mean()
    c = expected_payoff * np.exp(-r * T)     # in case r=0, this step is useless
    
    return c
#c = find_expected_payoff(s, k, r, ) # in case r=0, this step is useless
#c
df = pd.DataFrame(columns=['S0', 'K', 'V0', 'T', 'r', 'n', 'm', 'closed_solution', 'mc_price', 'ST_std'])
# fill the dataset
for s0 in S0:
    print(f's0 {s0}')
    for k in K:
        print(f'k {k}')
        for v0 in V0:
            print(f'v0 {v0}')
            for t in time_maturity:
                print(f't {t}')

                for r in int_rates:
                    # analytical closed sol
                    sol = call_price(kappa, theta, sigma, rho, v0, r, t, s0, k)

                    for i in range(len(n)):
                        nn = n[i]
                        mm = m[i]

                        # monte carlo solution:
                        ST = stock_price_generator(
                            t, nn, mm, r, s0, k, v0, sigma, theta, kappa, rho,
                            generator)  # stock prices of mc simulation

                        mc_price = find_expected_payoff(ST, k, r, t)

                        new_row = {
                            'S0': s0,
                            'K': k,
                            'V0': v0,
                            'T': t,
                            'V0': v0,
                            'r': r,
                            'n': nn,
                            'm': mm,
                            'closed_solution': sol,
                            'mc_price': mc_price,
                            'ST_std': np.std(ST)
                        }
                        
                        # append row to the dataframe
                        df = df.append(new_row, ignore_index=True)

df.head()
df.to_csv("solution_complete_mc.csv")
def plot_mc (s,T,n, m):
    # Setup figure
    plt.figure(figsize=(8,6))
    # noinspection PyTypeChecker
    ax_lines = plt.axes()

    maxx = np.max(s)
    minn = np.min(s)

    # Make the line plots
    t = np.linspace(0, T, num = n + 1)
    ns = 40
    for i in range(ns):
        ax_lines.plot(t, s[i, :], lw=1.0)
    ax_lines.set(xlabel='Years', ylabel='St',title='Price Simulations')
    ax_lines.set_xlim((0, T))
    ax_lines.set_ylim((minn, maxx))

    # Add mean value to line plots
    ax_lines.plot([0.0, 1.0], [s[:, -1].mean(), s[:, -1].mean()], lw='2', ls="--", label='mean')
    plt.show()

    # plot stock prices distributions
    plt.figure(figsize=(6,5))
    bins = np.arange(1.4, 3, .04)
    plt.hist(s[:, -1], bins=bins)
    plt.xlabel('ST')
    plt.ylabel('Number')
    plt.title('Stock Prices distribution')
    plt.show()
solution_mc = pd.read_csv('solution_complete_mc.csv')
solution_mc.head()
len(solution_mc)<jupyter_output><empty_output><jupyter_text>AV<jupyter_code>def stock_price_generator_av (T, n ,m, r, S0, k, V0, sigma, theta, kappa, rho, separate=False):
    dt = T / n
    
    m = int(m/2) # with av we can do half of the simulations
    n = int(n)
    
    # Brownian motions:
    dw_v = generator.normal(size=(m, n)) * np.sqrt(dt)
    dw_i = generator.normal(size=(m, n)) * np.sqrt(dt)

    dw_s = rho * dw_v + np.sqrt(1.0 - rho ** 2) * dw_i

    # Perform time evolution 
    s = np.empty((m, n + 1)) # initialisation stock prices vector
    s[:, 0] = S0
    
    s_ant = np.empty((m, n + 1))
    s_ant[:, 0] = S0

    v = np.ones(m) * V0
    v_ant = np.ones(m) * V0

    for t in range(n):
        dv = kappa * (theta - v) * dt + sigma * np.sqrt(v) * dw_v[:, t]
        dv_ant = kappa * (theta - v_ant) * dt + sigma * np.sqrt(v_ant) * dw_v[:, t]
        ds = r * s[:, t] * dt + np.sqrt(v) * s[:, t] * dw_s[:, t]
        ds_ant = r * s_ant[:, t] * dt + np.sqrt(v_ant) * s_ant[:, t] * dw_s[:, t]

        v = np.clip(v + dv, a_min=0.0, a_max=None)
        v_ant = np.clip(v_ant + dv_ant, a_min=0.0, a_max=None)
        
        s[:, t + 1] = s[:, t] + ds
        s_ant[:, t + 1] = s_ant[:, t] + ds_ant
        
    if separate==True:
        return s, s_ant
        
    return np.concatenate((s, s_ant))[:,-1]
mc_price_av = []
ST_av_std = []
i=0

for row in solution_mc.index:
    i = i+1
    if i % 10 == 0:
      print(i)

    T = solution_mc['T'][row]
    n = solution_mc['n'][row]
    m = solution_mc['m'][row]
    r = solution_mc['r'][row]
    S0 = solution_mc['S0'][row]
    K = solution_mc['K'][row]
    V0 = solution_mc['V0'][row]
    
    ST_ant = stock_price_generator_av (T, n, m, r, S0, K, V0, sigma, theta, kappa, rho)
    
    ST_av_std.append(np.std(ST_ant))
    mc_price_av.append(find_expected_payoff(ST_ant, K, r, T))

solution_mc['mc_price_av'] = mc_price_av
solution_mc['ST_av_std'] = ST_av_std
solution_mc.head()
standard_deviation = solution_mc[['ST_std', 'ST_av_std']]
standard_deviation.head()
prices = solution_mc[['closed_solution', 'mc_price', 'mc_price_av']]
prices.head()
rmse_mc = np.sqrt(((prices.closed_solution - prices.mc_price)**2).sum())
rmse_mc_av = np.sqrt(((prices.closed_solution - prices.mc_price_av)**2).sum())
rmse_mc
rmse_mc_av
solution_mc.to_csv("solution_complete_mc.csv")<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1586153846153846, 'fraction_numerical': 0.014769230769230769, 'mean_word_length': 2.042115114646701, 'pattern_counts': {'":': 0, '<': 17, '<?xml version=': 0, '>': 17, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '4', 'dir': 'jupyter-structured-clean-dedup', 'id': '62729', 'license': 'MIT', 'n_tokens_mistral': 2569, 'n_tokens_neox': 2493, 'n_words': 583, 'path': 'old (to ignore)/complete data/complete_mc.ipynb', 'repo_name': 'anitamezzetti/Heston_option_pricing', 'size': '6500'}