# Portfolio Optimization¶

In this tutorial, we consider portfolio optimization. A portfolio is a combination of financial instruments such as stocks and bonds. Unlike bank deposits, these financial instruments are risky assets with uncertain future payouts and returns. When investing in such risky assets, it is advisable to construct an appropriate portfolio and diversify investment across assets that are not highly correlated. In this way, if the value of one asset falls significantly, it can be covered by another asset that is not significantly correlated with its price.

In general, diversified investment involve a trade-off between risk and return. Portfolios must be constructed to properly reflect this trade-off, limiting return variability and risk of loss while maximizing return.

Here, we will use Fixstars Amplify to optimize the portfolio. The optimization uses estimates based on the historical data method. The historical data method is an estimation method that uses historical price data to determine assets' expected returns and risks. In this sample program, we will explain as follows

- Obtain historical data
- Formulate and implement a solver to find the optimal portfolio
- Run and evaluate the optimization
- Operation simulation

This sample program is intended to demonstrate an optimization application using Fixstars Amplify. The user is responsible for performing actual operation based on the information provided in this sample program.

## Obtain historical data¶

First, we implement the class `HistoricalDataLoader`

to get and use historical data.
`HistoricalDataLoader`

retrieves information of the necessary period from the stored data
based
on the given parameters, such as stocks, date of operation, operation period and sampling period.
This class uses `pandas_datareader`

to retrieve new historical data for the target stocks
from
the database in Stooq.

When executing this example program, it will take several minutes to download data from Stooq.

```
import datetime
import pathlib
import pandas as pd
from pandas_datareader import data
ONE_DAY = datetime.timedelta(days=1)
class HistoricalDataLoader:
def __init__(
self,
tickers: list[str], # List of tickers for stocks
operation_start_date: datetime.date, # Start date of operation
cache_dir: pathlib.Path, # Directory path of stored data
num_days_sampling=10, # Number of days of sampling period (10 business days for default)
num_days_operation=5, # Number of days of operation period (5 business days for default)
use_cache: bool = True, # Set `False` to download all data from Stooq
stdout: bool = False, # Set `True` for standard output
):
self._tickers = tickers
self._num_days_sampling = num_days_sampling
self._num_days_operation = num_days_operation
self._operation_start_date = operation_start_date
self._cache_dir = cache_dir
self._stdout = stdout
self._history = self._construct_history(use_cache)
self._sample_end_date: None | datetime.date = None
self._operation_end_date: None | datetime.date = None
def _construct_history(self, use_cache: bool) -> pd.DataFrame:
"""Construct historical data of all stocks"""
history = pd.DataFrame()
for ticker in self._tickers:
individual_history = self._fetch_save_individual_history(ticker, use_cache)
history = pd.concat([history, individual_history], axis=1).reindex(
index=individual_history.index
)
self._print_deletion_periods(history)
# Delete data on non-business days
return history.dropna(how="any")
def _print_deletion_periods(self, history: pd.DataFrame) -> None:
"""Print non-business days"""
mask = history.isnull().any(axis=1)
i = 1
while i < len(mask) - 1:
if mask.iloc[i] and not mask.iloc[i - 1]:
del_start_date = mask.index[i].to_pydatetime().date() # type: ignore
del_end_date = ""
for j in range(i, len(mask) - 1):
if not mask.iloc[j + 1]:
del_end_date = mask.index[j].to_pydatetime().date() # type: ignore
i = j
break
if self._stdout:
print(f"deletion {del_start_date}-{del_end_date}, nan is found.")
i += 1
def _calc_sample_end_date(self) -> datetime.date:
"""Calculate the end date of the sampling of historical data"""
check_date = self._operation_start_date
# The end date of the sampling is the business day closest to the start of operation before the start of operation.
# If you go back 14 days from the day before the start of the operation to the past, you should find at least one business day. (e.g., Chinese stocks around Chinese New Year)
for _ in range(14):
check_date -= ONE_DAY
key = check_date.strftime("%Y-%m-%d")
if key in self._history.index:
return check_date
raise RuntimeError(f"could not find a valid end date. {check_date=}")
def _calc_operation_end_date(self) -> datetime.date:
"""Calculate the end date of the operation period if the historical data include the operation period"""
hist = self.extract_operation_period()
if len(hist) == 0:
raise RuntimeError("historical data do not include operation period.")
return hist.index[0].to_pydatetime().date() # type: ignore
def _exist_operation_start_date(self) -> bool:
"""Return whether the historical data contain data on the start date of the operation"""
return self._operation_start_date.strftime("%Y-%m-%d") in self._history.index
def exist_operation_duration(self) -> bool:
"""Return whether the historical data contain data for the entire operational period"""
if not self._exist_operation_start_date():
return False
date_count = 0
date = self._operation_start_date
while date_count < self._num_days_operation:
if date.strftime("%Y-%m-%d") in self._history.index:
date_count += 1
date += ONE_DAY
if date > self._history.index[0].to_pydatetime().date(): # type: ignore
return False
return True
@property
def num_days_sampling(self) -> int:
"""Return number of days of sampling period"""
return self._num_days_sampling
@property
def num_days_operation(self) -> int:
"""Return number of days of operation period"""
return self._num_days_operation
@property
def tickers(self) -> list[str]:
"""Return list of tickers for stocks"""
return self._tickers
@property
def history(self) -> pd.DataFrame:
"""Return historical data for the entire period"""
return self._history
@property
def operation_end_date(self) -> datetime.date:
"""Return the end date of the operation period if the historical data include the operation period"""
if self._operation_end_date is None:
self._operation_end_date = self._calc_operation_end_date()
return self._operation_end_date
@property
def sample_end_date(self) -> datetime.date:
"""Return the end date of the sampling period"""
if self._sample_end_date is None:
self._sample_end_date = self._calc_sample_end_date()
return self._sample_end_date
def _load_cache(self, filepath: pathlib.Path, num_rows: int | None = None):
"""Load historical data from saved file"""
hist = pd.read_csv(filepath, index_col=0, nrows=num_rows)
hist.index = pd.to_datetime(hist.index) # type: ignore
if self._stdout:
print(f"Cache: {filepath.name} from {hist.index[-1]} to {hist.index[0]}") # type: ignore
return hist
def _download_history(self, ticker: str, start_date: datetime.date | None = None):
"""Download historical data from Stooq since `start_date`"""
if self._stdout:
print(f"Download: {ticker}, from {start_date}")
hist = data.DataReader(ticker, "stooq", start_date)
if len(hist) == 0:
return pd.DataFrame()
else:
return hist.drop(["Open", "High", "Low", "Volume"], axis=1)
def _save_cache(self, df: pd.DataFrame, filepath: pathlib.Path):
"""Save historical data to file"""
df.index.name = None
df.to_csv(filepath)
def _fetch_save_individual_history(
self, ticker: str, use_cache: bool = True
) -> pd.DataFrame:
"""Obtain and store historical data for individual stock"""
filepath = self._cache_dir.joinpath(ticker + ".csv")
start_date: datetime.date | None = None
hist = pd.DataFrame()
if filepath.exists() and use_cache:
hist = self._load_cache(filepath)
# Start date of the required download period
start_date = hist.index[0].date() + ONE_DAY # type: ignore
# +2 as a buffer to allow for a weekend
delta_days = datetime.timedelta(days=self._num_days_operation + 2)
if start_date is None or start_date <= self._operation_start_date + delta_days:
hist = pd.concat([self._download_history(ticker, start_date), hist])
self._save_cache(hist, filepath)
return hist.rename(columns={"Close": ticker}) # Consider only the closing price
def extract_sampling_period(self) -> pd.DataFrame:
"""Return historical data of the sampling period"""
i_end = self._history.index.get_loc(self.sample_end_date.strftime("%Y-%m-%d"))
ret = self._history.iloc[i_end : i_end + self._num_days_sampling] # type: ignore
return ret
def extract_operation_period(self) -> pd.DataFrame:
"""Return historical data of the operational period (if included in historical data)"""
key = self._operation_start_date.strftime("%Y-%m-%d")
i_start = self._history.index.get_loc(key) if key in self._history.index else 0
if not self.exist_operation_duration():
print("The historical data do not include the operation period.")
return pd.DataFrame()
ret = self._history.iloc[i_start - self._num_days_operation + 1 : i_start + 1] # type: ignore
return ret
```

## Formulate and implement a solver to find the optimal portfolio¶

We now formulate a formula to obtain the optimal portfolio using the estimates from the historical data method.

First, we assume the portfolio consists of $n$ financial instruments and will be operated for $d$ days. Let $w_i$ be the ratio of investment in a stock $i$ to the total investment, $w_i$ must satisfy the following equation:

$$ \sum_{i=1}^n w_i = 1 $$

Let $p_i$ be the closing price of a stock $i$, and the rate of return $r_i$ of a stock $i$ is defined using the value $p_{i,s}$ at the beginning of the investment and $p_{i,e}$ at the end of the investment as the following equation:

$$ r_i = (p_{i,e} - p_{i,s}) / p_{i,s} $$

In this case, the return on the entire portfolio $r_p$ is given by the following equation:

$$ r_p = \sum_{i=1}^{n} r_i w_i $$

On the other hand, the risk of diversified investment is usually the variance $\sigma^2_p$ or standard deviation $\sigma_p$ of the overall portfolio return, which can be described by the following equation:

$$ \sigma_p^2 = \sum_{i=1}^n \sum_{j=1}^n w_i w_j \sigma_{i,j}. $$

Here, $\sigma_{i,j}$ is the covariance of the returns of stocks $i$ and $j$. Generally, a portfolio is considered good if the return $R_p$ is large and the risk $\sigma_p^2$ is small. Therefore, we use the following mean and variance model $f(w_i)$ as the objective function, and optimization is performed to minimize it.

$$ f(w_i) = - r_p + \frac{\gamma}{2} \sigma_p^2. $$

Here, $\gamma$ is a parameter related to the balance between return and risk. In this sample program, $\gamma=1$ is the default setting for the portfolio optimization solver implemented below. If more emphasis is placed on risk reduction, adjust (increase) $\gamma$.

```
import amplify
import numpy as np
from typing import Any
import matplotlib.pyplot as plt
def calc_return_cov_rate(
h: HistoricalDataLoader, # Historical data
purchase_prices: np.ndarray, # Purchase prices of each stock
is_sampling_period: bool, # Boolean whether to use historical data of the sampling period
) -> tuple[np.ndarray | None, np.ndarray | None]:
"""Calculate the average return and covariance matrix for each stock over the operation period"""
if is_sampling_period:
prices_end = h.extract_sampling_period().to_numpy()
prices_start = np.roll(prices_end, -h.num_days_operation + 1, axis=0)
stop = h.num_days_sampling - h.num_days_operation + 1
else:
prices_end = h.extract_operation_period().to_numpy()
# If the data of the operation period is not included in the stored data
if len(prices_end) == 0:
return None, None
# The price of each stock at the start of the operation is the purchase price
prices_start = purchase_prices[np.newaxis, :]
# One return data set (i.e., no covariance matrix)
stop = 1
return_rates = ((prices_end - prices_start) / prices_start)[0:stop]
if len(return_rates) == 1 or len(return_rates[0]) == 1:
return return_rates.mean(axis=0), None
return return_rates.mean(axis=0), np.cov(return_rates, rowvar=False)
def calc_profit_risk(
h: HistoricalDataLoader, # Historical data
w_i_percent: amplify.PolyArray | np.ndarray, # Investment ratio (%)
purchase_prices: np.ndarray, # Purchase prices of stocks
is_sampling_period=True, # Boolean whether to use historical data of the sampling period
) -> tuple[Any, Any]:
"""Calculate and return the return and variance for the entire given portfolio"""
w_i = w_i_percent / 100
return_rates, cov_rates = calc_return_cov_rate(
h, purchase_prices, is_sampling_period
)
if return_rates is None:
return None, None
profit = amplify.einsum("i,i->", w_i, return_rates)
if cov_rates is None:
return profit, None
risk = amplify.einsum("i,j,ij->", w_i, w_i, cov_rates)
return profit, risk
def generate_objective(
h: HistoricalDataLoader, # Historical data
w_i_percent: amplify.PolyArray, # Investment ratio (%)
purchase_prices: np.ndarray, # Purchase prices of stocks
gamma: float, # Gamma in mean and variance model
) -> amplify.Poly:
"""Construct and return an objective function based on the mean/variance model and an additive term"""
profit, risk = calc_profit_risk(h, w_i_percent, purchase_prices)
return -profit + gamma * 0.5 * risk
def generate_constraint(w_i_percent: amplify.PolyArray) -> amplify.Constraint:
"""Construct and return the constraint that results in a total investment ratio of 100%"""
return amplify.equal_to(w_i_percent, 100)
def solve(
h: HistoricalDataLoader, # Historical data
purchase_prices: np.ndarray, # Purchase prices of stocks
gamma=1, # Gamma in mean and variance model
constraint_weight: float = 1, # Weight for constraint condition
timeout_ms: int = 5000, # Solving timeout (ms)
num_solves: int = 3, # Number of serial runs for solve (https://amplify.fixstars.com/ja/docs/amplify/v1/serial.html)
plot_optim_history=False,
) -> np.ndarray:
"""Optimize portfolio based on historical data, estimated purchase price of each stock, and investment budget"""
gen = amplify.VariableGenerator()
# q: Investment ratio (%)
# For diversification, limit ratio of each stock to 50%
q = gen.array("Integer", len(h.tickers), bounds=(0, 50))
objective = generate_objective(h, q, purchase_prices, gamma)
constraint = generate_constraint(q)
model = amplify.Model(objective, constraint_weight * constraint)
client = amplify.FixstarsClient()
client.parameters.timeout = datetime.timedelta(milliseconds=timeout_ms)
client.parameters.outputs.num_outputs = 0 # Return all solutions found
# client.token = "xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx" # If you use Amplify in a local environment or Google Colaboratory, enter your Amplify API token.
# Solve the problem
result = amplify.solve(model, client, num_solves=num_solves)
if len(result) == 0:
print("no solution is found")
if plot_optim_history:
optim_history(result)
return q.evaluate(result.best.values)
def optim_history(result: amplify.Result) -> None:
"""Display optimization history"""
ax = plt.figure().add_subplot()
# Obtain the time of each solution and the value of the objective function
for i in range(result.num_solves):
times = [solution.time.total_seconds() for solution in result.split[i]]
objectives = [solution.objective for solution in result.split[i]]
ax.plot(times, objectives, marker="*")
ax.set_xlabel("elapsed time in seconds")
ax.set_ylabel("objective value")
ax.grid(True)
plt.show()
```

## Run and evaluate the optimization¶

We will now construct the optimal portfolio for the stocks that make up the NASDAQ 100 and then compare the constructed optimal portfolio to the profit margin of the NASDAQ 100 Index. Here, the NASDAQ 100 Index is simply a portfolio based on the market capitalization of each component.

As a problem set, the historical data used for portfolio optimization is ten business days back from
the
start of the operation, and the operation period is five business days.
These values are set as defaults in the constructor of the `HistoricalDataLoader`

class.

In this sample program, January 1, 2024 is assumed as the start date of operation. This is for verification after the operation period. In typical portfolio optimization, the start date is the date on which you begin working with the optimized portfolio.

Here, historical data for dates that are not business days are deleted at the time of the
`HistoricalDataLoader`

class instantiation.
A business day is defined as "a day on which all stocks under consideration can be traded."

### Obtain historical data¶

Now, let's use the `HistoricalDataLoader`

to get historical data for the stocks that make
up
the NASDAQ100.

```
# Directory path for storing historical data
history_data_dir = pathlib.Path("../../../storage/portfolio")
# Start date of operation
date_of_operation = datetime.date.fromisoformat("2024-04-01")
# NASDAQ100 component tickers (May 18, 2024)
nasdaq100_stocks = ["ADBE", "ADP", "ABNB", "GOOGL", "GOOG", "AMZN", "AMD", "AEP", "AMGN", "ADI",
"ANSS", "AAPL", "AMAT", "ASML", "AZN", "TEAM", "ADSK", "BKR", "BIIB", "BKNG",
"AVGO", "CDNS", "CDW", "CHTR", "CTAS", "CSCO", "CCEP", "CTSH", "CMCSA", "CEG",
"CPRT", "CSGP", "COST", "CRWD", "CSX", "DDOG", "DXCM", "FANG", "DLTR", "DASH",
"EA", "EXC", "FAST", "FTNT", "GEHC", "GILD", "GFS", "HON", "IDXX", "ILMN",
"INTC", "INTU", "ISRG", "KDP", "KLAC", "KHC", "LRCX", "LIN", "LULU", "MAR",
"MRVL", "MELI", "META", "MCHP", "MU", "MSFT", "MRNA", "MDLZ", "MDB", "MNST",
"NFLX", "NVDA", "NXPI", "ORLY", "ODFL", "ON", "PCAR", "PANW", "PAYX", "PYPL",
"PDD", "PEP", "QCOM", "REGN", "ROP", "ROST", "SIRI", "SBUX", "SNPS", "TTWO",
"TMUS", "TSLA", "TXN", "TTD", "VRSK", "VRTX", "WBA", "WBD", "WDAY", "XEL", "ZS"] # fmt: skip
# Instantiate historical data class
nasdaq100_stocks_history = HistoricalDataLoader(
nasdaq100_stocks,
date_of_operation,
history_data_dir,
stdout=True,
)
```

As a test, let's use `nasdaq100_stocks_history`

to display historical data for Adobe
Corporation (`"ADBE"`

).
Data is displayed for the ten business days back to the business day before the start of operation.

```
nasdaq100_stocks_history.extract_sampling_period()["ADBE"]
```

### Determine the purchase prices¶

During market hours, the value of stocks changes from moment to moment. The purchase price of a stock, whether limit price or market order, is not known until it is purchased. However, the purchase price must be considered to determine which and how many stocks to buy in an optimal portfolio and evaluate that portfolio. We assume that each stock can be purchased at a randomly determined cost between 0 and 1% more than the previous day's closing market price on the day the portfolio is launched. The premiums are also assumed to take into account commissions and other fees.

```
rng = np.random.default_rng(0)
def calc_purchase_prices(h: HistoricalDataLoader, low=1.0, high=1.01) -> np.ndarray:
"""Calculate the purchase prices of the subject stocks. Add random prices in the range of [low, high] to each closing price on the day before the operation start date."""
# Closing prices on the day before the operation start date (i.e., the last day of the sampling period)
close = h.extract_sampling_period()[h.sample_end_date : h.sample_end_date].iloc[0]
return close.to_numpy() * rng.uniform(low, high, size=len(close))
# The purchase prices at start of operation
purchase_prices_stocks = calc_purchase_prices(nasdaq100_stocks_history)
```

### Solve for optimal portfolios¶

Now, we will find the optimal portfolio based on the historical data with the solver implemented above.

```
solution = solve(nasdaq100_stocks_history, purchase_prices_stocks)
```

### Evaluate the portfolio¶

Now, we evaluate the obtained portfolio.
In the following `evaluate()`

, (1) the rate of return and risk for the historical data, and
(2)
the rate of return after the investment (if data during the investment period is given) are calculated
for
the portfolio to be evaluated.
After subtracting the tax on the gain (20% separate tax in Japan), the final rate of return is
returned as
the return value.

```
import math
TAX_RATE = 0.2
def evaluate(
sol: np.ndarray,
purchase_prices: np.ndarray,
h: HistoricalDataLoader,
show_historical=False,
tax_rate=TAX_RATE,
) -> float | None:
"""Evaluate the investment results based on the portfolio and return the overall portfolio rate of return based on historical data and (if data exist) the rate of return over the investment period."""
portfolio = {
ticker: f"{int(sol[i])}%" for i, ticker in enumerate(h.tickers) if sol[i] > 0
}
print(f"portfolio: {portfolio}")
# Rate of return/risk for the sampling period
h_profit, h_risk = calc_profit_risk(h, sol, purchase_prices)
if show_historical:
print(f"historical profit: {float(h_profit)*100:.1f} %")
if h_risk is not None:
# Display standard deviation value as risk
print(f"historical risk: {math.sqrt(float(h_risk))*100:.1f} %")
op_profit, _ = calc_profit_risk(h, sol, purchase_prices, is_sampling_period=False)
if op_profit is None:
return None
# Rate of return after operation (Only if data during the investment period is given, `tax_rate` on earnings is also considered)
print(f"operational profit: {float(op_profit)*100:.1f} %")
return (
float(op_profit) * (1 - tax_rate) if float(op_profit) > 0 else float(op_profit)
)
# Evaluate the optimal portfolio found in the upper cell
return_rate = evaluate(
solution,
purchase_prices_stocks,
nasdaq100_stocks_history,
show_historical=True,
)
```

Next, we will also evaluate the NASDAQ 100 Index using the same starting date and investment period as the portfolio construction conditions. In the NASDAQ 100 Index, constituent stocks are structured proportionately to market capitalization.

```
nasdaq100_index_history = HistoricalDataLoader(
["^NDX"], # Ticker of NASDAQ 100 Index
date_of_operation,
history_data_dir,
)
# The purchase price at start of operation
purchase_prices_index = calc_purchase_prices(nasdaq100_index_history)
# The investment ratio is set to 100% since only one stock (index) is used
return_rate = evaluate(
np.array([100]),
purchase_prices_index,
nasdaq100_index_history,
show_historical=True,
)
```

## Operation simulation using portfolio optimization solver¶

Now, we will perform an operation simulation using the portfolio optimization solver implemented in
this
sample program.
The first business day after January 1, 2024, is set as the operation start date, and from this date,
five
operations are performed for each sampling period of historical data.
The operation and optimization conditions are the same as the default values of
`HistoricalDataLoader`

, with the operation period as `num_days_operation = 5`

business days, and `num_days_sampling = 10`

business days back from each operation date as
the
historical data to be used for optimization on each operation date.
Based on the portfolio to be constructed in each investment round, the funds increased or decreased by
the
previous investment are invested, including compound interest.
The following figure is a schematic of such an investment flow.

```
def simulate_operation(
tickers: list[str], # Tickers of stocks
operation_date: datetime.date | list[datetime.date], # Start date of operation
num_operation_rounds: int, # Number of operation rounds
) -> tuple[list[datetime.date], list[float]]:
"""Operation simulation of buying and selling over `num_operation_rounds` time's investment rounds."""
# Scan historical data and store the start and end dates of each operation round in `operation_dates`
operation_dates: list[datetime.date] = []
if isinstance(operation_date, list):
operation_dates = operation_date
else:
i_operation_round = 0
while i_operation_round < num_operation_rounds:
while True:
h = HistoricalDataLoader(tickers, operation_date, history_data_dir)
# If `operation_date` is not a business day, revalidate one day later as a business day
if not h.exist_operation_duration():
operation_date += ONE_DAY
break
i_operation_round += 1
operation_dates.append(operation_date)
operation_dates.append(h.operation_end_date)
operation_date += datetime.timedelta(days=h.num_days_sampling)
break
# List to include the amount of assets at the beginning and end of each investment round
fund_evolution: list[float] = []
funds = 1.0 # Initial asset value
print(f"Initial funds: {funds:.2f}")
for i, date in enumerate(operation_dates):
# If `i` is odd, `date` is the end date of the operation
if i % 2 == 1:
continue
h = HistoricalDataLoader(tickers, date, history_data_dir)
print("============================================")
print(f"Operation round: {int(i / 2) + 1} of {num_operation_rounds}, {date}")
# Add the amount of assets at the start of the current round to `fund_evolution`
fund_evolution.append(funds)
# Purchase prices of all stocks at the start of current round
prices = calc_purchase_prices(h)
# If there is more than one candidate of stocks, find the optimal portfolio; if one candidate, buy 100% of that stock
if len(prices) > 1:
solution = solve(h, prices)
else:
solution = np.array([100])
# Evaluate the investment results in current round and calculate the amount of assets after increase/decrease
operational_return = evaluate(solution, prices, h)
if operational_return is not None:
funds *= 1.0 + operational_return
print(f"Current funds: {funds:.2f}")
# Add the amount of assets at the end of the current round to `fund_evolution`
fund_evolution.append(funds)
return operation_dates, fund_evolution
```

```
start_date = datetime.date.fromisoformat(
"2024-01-01"
) # Start date of operation of the first round
num_operation_rounds = 5 # Number of operation rounds
operation_dates, optim_portfo_evolution = simulate_operation(
nasdaq100_stocks, start_date, num_operation_rounds
)
```

The following cell is a simulation of an operation in which stocks based on the NASDAQ 100 Index, rather than the optimal portfolio, are traded on the same day as the above investment round.

```
_, index_op_evolution = simulate_operation(["^NDX"], start_date, num_operation_rounds)
```

The simulation below shows a case where stocks are purchased based on the NASDAQ100 Index on the start date of the first investment round and are not traded until the end date of the last investment round.

```
round_start_date = operation_dates[0]
round_end_date = operation_dates[-1]
h = HistoricalDataLoader(["^NDX"], round_start_date, history_data_dir)
# Purchase price on the start date of the first investment round
prices = calc_purchase_prices(h)
# Index on the last day of the last operational round
i_end = h.history.index.get_loc(round_end_date.strftime("%Y-%m-%d"))
# Taxed funds after sale on the last day of the last investment round
return_rate = h.history.iloc[i_end].iloc[0] / prices[0]
if return_rate > 1:
return_rate = (return_rate - 1) * (1 - TAX_RATE) + 1
index_no_op_evolution = [1, return_rate]
```

Finally, we plot the results of operating the components of the NASDAQ 100 Index in the following three ways.

`optim portfo`

: Operate according to the optimal portfolio`index-op`

: Operate according to the index`index-no-op`

: Buy based on the index and do not buy or sell until the end

It is important to note that since the stocks are purchased at a random premium of 0~1%,
`optim portfo`

and `index-op`

, which repeat the buying and selling for each
investment round, will be affected negatively by the surcharges and taxes.
The aggressive trading of financial instruments based on the optimal portfolio is particularly
detrimental
when the market is downtrend.

In this sample program, the sampling period and the operation period of the historical data are set
to 10
days and 5 days, respectively, as default values of the constructor of the historical data class
`HistoricalDataLoader`

.
If the price fluctuation over a period of time is relatively small for the portfolio as a whole, it is
possible to reduce the number of trades by extending the operation period.

In addition, a lower-risk, higher-return portfolio can be expected by adding a more comprehensive
range
of financial instruments other than stocks (e.g., metals, fossil fuels, cryptocurrencies) to the
portfolio.
For example, in Stooq used in this sample program, gold can be obtained with the ticker
`GC.C`

,
crude oil with `CL.C`

, and bitcoin/USD with `BTCUSD`

.
Let's add various commodities to your portfolio and run the investment simulation.

```
import matplotlib.pyplot as plt
from matplotlib import dates as mdates
ax = plt.figure().add_subplot()
color_stocks = [0.2, 0.2, 1]
color_index = [1, 0.2, 0.2]
# Plot "optim portfo", "index-op" for all operational rounds
operation_dates = np.array(operation_dates)
ax.plot(operation_dates, optim_portfo_evolution, linestyle=":", color=color_stocks)
ax.plot(operation_dates, index_op_evolution, linestyle=":", color=color_index)
for i, date in enumerate(operation_dates):
text = "E"
if i % 2 == 0:
text = "S"
(optimal,) = ax.plot(
operation_dates[i : i + 2],
optim_portfo_evolution[i : i + 2],
color=color_stocks,
marker="*",
)
(index,) = ax.plot(
operation_dates[i : i + 2],
index_op_evolution[i : i + 2],
color=color_index,
marker="*",
)
ax.annotate(text, (operation_dates[i], optim_portfo_evolution[i]), fontsize=10)
ax.annotate(text, (operation_dates[i], index_op_evolution[i]), fontsize=10)
# Plot "index-no-op" (Only on the first day of the first operation round and the last day of the last operation round)
(no_transact,) = ax.plot(
[operation_dates[0], operation_dates[-1]],
index_no_op_evolution,
marker="o",
color=color_index,
)
ax.legend(
[optimal, index, no_transact],
["optim portfo", "index-op", "index-no-op"],
loc="lower right",
)
ax.annotate(
" S: Start date of round\n E: End date",
(operation_dates[0], max(optim_portfo_evolution) * 0.99),
fontsize=10,
)
ax.set_xlabel("Date", fontsize=10)
ax.set_ylabel("Total asset", fontsize=10)
ax.tick_params(labelsize=10)
ax.set_xlim((operation_dates[0] - ONE_DAY, operation_dates[-1] + ONE_DAY))
ax.xaxis.set_major_formatter(mdates.DateFormatter("%m/%d"))
plt.show()
```