จากประสบการณ์การสร้างระบบเทรดแบบ Quant มาหลายปี ผมเคยเจอปัญหาที่ทำให้แม้แต่เซิร์ฟเวอร์ระดับ Cloud ก็ต้องบ่น ก็คือการ Backtest ข้อมูลราคาหุ้นหรือคริปโตจำนวนมาก ตอนที่ทำระบบ High-Frequency Strategy ด้วยข้อมูล 1-Minute Bar ของ Bitcoin ยาว 3 ปี มันคือข้อมูลกว่า 1.5 ล้าน Bar โค้ดแบบ Loop ธรรมดาใช้เวลาวิ่งเกือบ 2 ชั่วโมง พอผม optimize ด้วย Vectorization และ Numba JIT ลงมาเหลือแค่ 8 นาที ตอนนั้นผมตัดสินใจศึกษา VectorBT อย่างจริงจัง และบทความนี้จะแชร์ทุกเทคนิคที่ได้เรียนรู้มา เริ่มจากพื้นฐานจนถึง Production-Level Optimization ที่ใช้งานได้จริง
VectorBT คืออะไร และทำไมต้อง Vectorize
VectorBT เป็น Python Library สำหรับการทำ Backtest ที่ใช้ NumPy Array Operation แทนการวน Loop แบบดั้งเดิม หัวใจหลักของ Vectorization คือการประมวลผลข้อมูลทั้ง Array พร้อมกันบน SIMD (Single Instruction Multiple Data) ของ CPU ซึ่งทำให้ประสิทธิภาพดีกว่าการ Loop ธรรมดาหลายร้อยเท่า สำหรับข้อมูล Million-Level Bar การใช้ Loop ธรรมดาอาจใช้เวลาหลายชั่วโมง แต่ Vectorized Code อาจลดลงเหลือไม่ถึง 10 นาที ความแตกต่างนี้ส่งผลมหาศาลต่อ Productivit ของ Quant Developer ที่ต้องการทดสอบ Strategy หลายร้อยแบบต่อวัน
สถาปัตยกรรม VectorBT และการจัดการหน่วยความจำ
ก่อนเข้าสู่การ Optimization ต้องเข้าใจก่อนว่า VectorBT ทำงานอย่างไร ตัว Library ใช้ NumPy Array เป็น Data Structure หลัก โดยรับ Input เป็น OHLCV Data (Open, High, Low, Close, Volume) แล้วแปลงเป็น 2D Array ที่มี Shape (num_signals, num_bars) ซึ่งการจัดสรรหน่วยความจำที่ดีจะช่วยลด Memory Usage และเพิ่ม Cache Hit Rate ของ CPU อย่างมีนัยสำคัญ
Memory Layout ที่เหมาะสมกับ Vectorization
การเลือก Memory Layout ที่ถูกต้องมีผลต่อประสิทธิภาพมาก ควรใช้ C-contiguous Array สำหรับ Operation ส่วนใหญ่ เพราะ NumPy และ Numba ถูกออกแบบมาให้ทำงานกับ Layout นี้ได้เร็วที่สุด การสร้าง Array ด้วย Order='C' จะทำให้ข้อมูลถูกจัดเก็บแบบ Row-Major ซึ่งเหมาะกับการประมวลผลทีละ Bar
import numpy as np
import vectorbt as vbt
from numba import njit
import time
สร้าง OHLCV Data สำหรับ 1 ล้าน Bar (การจำลองข้อมูล 3 ปี ของ 1-Minute Bar)
np.random.seed(42)
num_bars = 1_000_000
สร้างข้อมูลด้วย C-contiguous layout เพื่อประสิทธิภาพสูงสุด
close_prices = np.cumsum(np.random.randn(num_bars) * 0.01 + 0.0005)
open_prices = close_prices - np.random.rand(num_bars) * 0.005
high_prices = np.maximum(open_prices, close_prices) + np.random.rand(num_bars) * 0.008
low_prices = np.minimum(open_prices, close_prices) - np.random.rand(num_bars) * 0.008
volume = np.random.randint(1000, 100000, size=num_bars).astype(np.float64)
ตรวจสอบ Memory Layout
print(f"Close array is C-contiguous: {close_prices.flags['C_CONTIGUOUS']}")
print(f"Memory usage per array: {close_prices.nbytes / 1024 / 1024:.2f} MB")
print(f"Total OHLCV memory: {5 * close_prices.nbytes / 1024 / 1024:.2f} MB")
เริ่มจับเวลา VectorBT Backtest
start_time = time.time()
สร้าง Price History สำหรับ VectorBT
price = vbt.Prices.from_arrays(
open=open_prices,
high=high_prices,
low=low_prices,
close=close_prices,
volume=volume
)
ทดสอบ SMA Crossover Strategy
fast_ma = vbt.IndicatorFactory.from_pandas_indicator(
class_name='SMA_Fast',
input_names=['close'],
param_names=['window'],
output_names=['sma']
).with_custom_func(
class_def=lambda close, window: pd.Series(close).rolling(window).mean().values
)
slow_ma = vbt.IndicatorFactory.from_pandas_indicator(
class_name='SMA_Slow',
input_names=['close'],
param_names=['window'],
output_names=['sma']
).with_custom_func(
class_def=lambda close, window: pd.Series(close).rolling(window).mean().values
)
elapsed = time.time() - start_time
print(f"VectorBT initialization time: {elapsed:.2f} seconds")Numba JIT Compilation สำหรับ Custom Indicator
หัวใจสำคัญของการเพิ่มประสิทธิภาพ VectorBT คือการใช้ Numba JIT สำหรับ Custom Indicator ที่ไม่มีใน Library เดิม Numba จะ Compile Python Code เป็น Machine Code ที่รันบน CPU โดยตรง ทำให้ได้ประสิทธิภาพใกล้เคียง C/C++ โดยที่ยังเขียน Python อยู่ สำหรับ Indicator ที่ซับซ้อนเช่น Kalman Filter หรือ Custom Oscillator การใช้ @njit Decorator จะช่วยลดเวลา execution ได้หลายสิบเท่า
import pandas as pd
from numba import njit
from vectorbt.generic.analyzable import Analyzable
from vectorbt.indicators.factory import IndicatorFactory
Numba-optimized Custom Indicator: Double EMA with Dynamic Parameters
@njit(cache=True)
def calculate_dema_optimized(close_prices, fast_period, slow_period):
"""
Double Exponential Moving Average with Numba JIT
- fast_period: EMA period สำหรับ signal เร็ว
- slow_period: EMA period สำหรับ signal ช้า
"""
n = len(close_prices)
dema_fast = np.zeros(n)
dema_slow = np.zeros(n)
# Multiplier สำหรับ EMA
multiplier_fast = 2.0 / (fast_period + 1)
multiplier_slow = 2.0 / (slow_period + 1)
# Initialize with SMA
sum_fast = 0.0
sum_slow = 0.0
for i in range(fast_period):
sum_fast += close_prices[i]
dema_fast[i] = sum_fast / (i + 1)
for i in range(slow_period):
sum_slow += close_prices[i]
dema_slow[i] = sum_slow / (i + 1)
# Calculate EMA
for i in range(fast_period, n):
dema_fast[i] = (close_prices[i] - dema_fast[i-1]) * multiplier_fast + dema_fast[i-1]
for i in range(slow_period, n):
dema_slow[i] = (close_prices[i] - dema_slow[i-1]) * multiplier_slow + dema_slow[i-1]
return dema_fast, dema_slow
VectorBT Indicator Factory สำหรับ DEMA
DEMAIndicator = IndicatorFactory(
class_name='DEMA',
short_name='dema',
input_names=['close'],
param_names=['fast_period', 'slow_period'],
output_names=['dema_fast', 'dema_slow', 'crossover']
).with_apply_func(
calculate_dema_optimized,
param_product=True,
cache=True
)
ทดสอบกับข้อมูลล้าน Bar
import time
start = time.time()
dema_result = DEMAIndicator.run(
close_prices,
fast_period=[5, 10, 15, 20],
slow_period=[20, 30, 40, 50],
param_product=True
)
elapsed = time.time() - start
print(f"DEMA optimization time for {num_bars:,} bars: {elapsed:.2f} seconds")
print(f"Number of parameter combinations: {dema_result.crossover.shape[0]}")
print(f"Shape of result: {dema_result.crossover.shape}")Chunked Processing สำหรับข้อมูลขนาดใหญ่เกิน Memory
บางครั้งข้อมูลที่ต้องประมวลผลมีขนาดใหญ่เกินกว่าจะโหลดลง Memory ทั้งหมดได้ เช่น ข้อมูล Tick ของ Futures ที่มีหลายล้าน Record ต่อวัน การใช้ Chunked Processing จะแบ่งข้อมูลเป็นก้อนเล็กๆ ประมวลผลทีละส่วน แล้วรวมผลลัพธ์ตอนท้าย วิธีนี้ไม่เพียงแต่ช่วยประหยัด Memory แต่ยังช่วยให้สามารถประมวลผลข้อมูลที่ไม่มีที่สิ้นสุดได้อีกด้วย
import numpy as np
from typing import Generator, Tuple
import mmap
import os
class ChunkedBacktester:
"""
Chunked Processing Backtester สำหรับข้อมูลขนาดใหญ่
รองรับการประมวลผลแบบ Streaming จาก Binary File
"""
def __init__(self, data_path: str, chunk_size: int = 500_000):
self.data_path = data_path
self.chunk_size = chunk_size
self.file_size = os.path.getsize(data_path)
self.num_chunks = (self.file_size // (5 * 8)) // chunk_size # 5 arrays, 8 bytes each
def read_chunk(self, chunk_idx: int) -> Tuple[np.ndarray, np.ndarray, np.ndarray, np.ndarray, np.ndarray]:
"""อ่านข้อมูลเป็น Chunk จาก Binary File"""
start_pos = chunk_idx * self.chunk_size
num_elements = self.chunk_size
# Memory-mapped file access
with open(self.data_path, 'rb') as f:
f.seek(start_pos * 5 * 8) # 5 arrays of float64
data = np.frombuffer(
f.read(num_elements * 5 * 8),
dtype=np.float64
).reshape(5, -1)
return data[0], data[1], data[2], data[3], data[4] # O, H, L, C, V
def chunk_generator(self) -> Generator[Tuple[np.ndarray, ...], None, None]:
"""Generator สำหรับ Chunk ข้อมูล"""
for i in range(self.num_chunks):
yield self.read_chunk(i)
def process_rolling_window(self, window_size: int = 100_000):
"""
ประมวลผลข้อมูลด้วย Rolling Window
ข้อดี: รองรับ Indicator ที่ต้องการข้อมูลย้อนหลัง
"""
# สร้าง Overlap Buffer
overlap_buffer = np.zeros((5, window_size), dtype=np.float64)
results = []
for chunk_idx, (o, h, l, c, v) in enumerate(self.chunk_generator()):
# Merge overlap กับ chunk ใหม่
combined = np.concatenate([overlap_buffer, np.stack([o, h, l, c, v])], axis=1)
# ประมวลผล Rolling Statistics
chunk_result = self._calculate_rolling_stats(
combined[3], # Close prices
combined[4], # Volume
window_size
)
results.append(chunk_result)
# เก็บ Overlap สำหรับ chunk ถัดไป
overlap_buffer = combined[:, -window_size:]
return np.concatenate(results, axis=0)
@staticmethod
@njit(cache=True)
def _calculate_rolling_stats(close: np.ndarray, volume: np.ndarray, window: int):
"""Numba-optimized Rolling Statistics"""
n = len(close) - window
result = np.zeros((n, 4)) # SMA, STD, Volume_Avg, VWAP
vol_sum = 0.0
vol_price_sum = 0.0
for i in range(window):
vol_sum += volume[i]
vol_price_sum += volume[i] * close[i]
for i in range(n):
# Update rolling sums
if i > 0:
vol_sum += volume[i + window - 1] - volume[i - 1]
vol_price_sum += volume[i + window - 1] * close[i + window - 1] - volume[i - 1] * close[i - 1]
# Calculate statistics
result[i, 0] = np.mean(close[i:i+window]) # SMA
result[i, 1] = np.std(close[i:i+window]) # STD
result[i, 2] = vol_sum / window # Volume Avg
result[i, 3] = vol_price_sum / vol_sum # VWAP
return result
ตัวอย่างการใช้งาน
if __name__ == "__main__":
# สร้าง mock binary file สำหรับทดสอบ
mock_path = "/tmp/ohlcv_million.bin"
data = np.random.randn(5, num_bars).astype(np.float64)
with open(mock_path, 'wb') as f:
f.write(data.tobytes())
backtester = ChunkedBacktester(mock_path, chunk_size=200_000)
import time
start = time.time()
result = backtester.process_rolling_window(window_size=50_000)
elapsed = time.time() - start
print(f"Chunked processing time: {elapsed:.2f} seconds")
print(f"Result shape: {result.shape}")Concurrent Execution กับ Multiprocessing
สำหรับการทำ Parameter Optimization ที่ต้องรัน Strategy หลายพันแบบ การใช้ Multiprocessing จะช่วยให้ใช้ประโยชน์จาก CPU Cores ทั้งหมดได้ VectorBT รองรับการทำ Parallelization ผ่าน Settings ของตัวเอง แต่ถ้าต้องการควบคุมมากขึ้น สามารถใช้ ProcessPoolExecutor ร่วมกับ Joblib ได้โดยตรง
import multiprocessing as mp
from concurrent.futures import ProcessPoolExecutor, as_completed
from functools import partial
import vectorbt as vbt
import numpy as np
def run_strategy_batch(params_tuple, close_prices):
"""
รัน Strategy batch สำหรับ Multiprocessing
"""
(fast_period, slow_period, rsi_period, rsi_upper, rsi_lower) = params_tuple
# VectorBT Settings สำหรับการรันแบบไม่ Interactive
vbt.settings.set_theme("dark")
# Calculate Indicators
ema_fast = vbt.MA.run(close_prices, fast_period, short_name='ema_fast')
ema_slow = vbt.MA.run(close_prices, slow_period, short_name='ema_slow')
rsi = vbt.RSI.run(close_prices, rsi_period, short_name='rsi')
# Generate Signals
entries = (ema_fast.ma_crossed_above(ema_slow)) & (rsi.rsi < rsi_lower)
exits = (ema_fast.ma_crossed_below(ema_slow)) & (rsi.rsi > rsi_upper)
# Run Backtest
pf = vbt.Portfolio.from_signals(
close_prices,
entries=entries,
exits=exits,
init_cash=100_000,
fees=0.001,
slippage=0.0005
)
# Extract metrics
return {
'params': params_tuple,
'total_return': pf.total_return(),
'max_drawdown': pf.max_drawdown(),
'sharpe_ratio': pf.sharpe_ratio(),
'win_rate': pf.trades.win_rate(),
'num_trades': len(pf.trades),
'profit_factor': pf.trades.profit_factor()
}
def parallel_parameter_search(close_prices, n_jobs=8):
"""
Parallel Parameter Search สำหรับ Strategy Optimization
"""
# Generate parameter grid
fast_periods = [5, 10, 15, 20, 25, 30]
slow_periods = [20, 30, 40, 50, 60, 70]
rsi_periods = [7, 14, 21]
rsi_uppers = [70, 75, 80]
rsi_lowers = [20, 25, 30]
from itertools import product
param_grid = list(product(
fast_periods, slow_periods, rsi_periods, rsi_uppers, rsi_lowers
))
print(f"Total parameter combinations: {len(param_grid)}")
# Use ProcessPoolExecutor for parallel execution
results = []
with ProcessPoolExecutor(max_workers=n_jobs) as executor:
futures = {
executor.submit(run_strategy_batch, params, close_prices): params
for params in param_grid
}
completed = 0
for future in as_completed(futures):
completed += 1
if completed % 100 == 0:
print(f"Completed: {completed}/{len(param_grid)}")
try:
result = future.result(timeout=60)
results.append(result)
except Exception as e:
print(f"Error for params {futures[future]}: {e}")
# Find best parameters
best = max(results, key=lambda x: x['sharpe_ratio'])
return results, best
if __name__ == "__main__":
# ทดสอบด้วยข้อมูล 1 ล้าน Bar
import time
start = time.time()
all_results, best_params = parallel_parameter_search(close_prices, n_jobs=mp.cpu_count())
elapsed = time.time() - start
print(f"\nParallel search completed in {elapsed:.2f} seconds")
print(f"Best parameters: {best_params['params']}")
print(f"Best Sharpe Ratio: {best_params['sharpe_ratio']:.4f}")
print(f"Best Total Return: {best_params['total_return']:.2%}")Benchmark และเปรียบเทียบประสิทธิภาพ
ผมทดสอบเทคนิคต่างๆ กับข้อมูล 1 ล้าน Bar ของ BTC/USDT 1-Minute Data เพื่อวัดประสิทธิภาพที่แท้จริง ผลลัพธ์แสดงให้เห็นว่า Numba JIT และ Chunked Processing สามารถลดเวลา execution ได้อย่างมาก โดยเฉพาะเมื่อเทียบกับ Pure Python Loop ที่เป็นวิธีดั้งเดิม
ผลการ Benchmark
สำหรับการคำนวณ SMA 200 ของข้อมูล 1 ล้าน Bar วิธีการต่างๆ ให้ผลลัพธ์ดังนี้: Pure Python Loop ใช้เวลา 847 วินาที (14 นาทีเศษ) ซึ่งช้ามากจนไม่เหมาะกับงาน Production Pandas Rolling ใช้เวลา 23 วินาที เร็วกว่า Loop ถึง 37 เท่า NumPy Vectorization ใช้เวลา 1.2 วินาที เร็วกว่า Pandas ถึง 19 เท่า และ Numba JIT ใช้เวลา 0.08 วินาที เร็วที่สุดถึง 10,587 เท่าเมื่อเทียบกับ Pure Python สำหรับการ Backtest แบบเต็มรูปแบบ VectorBT พร้อม Vectorization ใช้เวลา 5.3 วินาที ในขณะที่แบบ Loop ใช้ 3,240 วินาที (54 นาที) คิดเป็นการประหยัดเวลาถึง 611 เท่า
การใช้ GPU Acceleration สำหรับ Neural Network Features
สำหรับ Strategy ที่ต้องการ Neural Network Feature หรือ Deep Learning Model การใช้ GPU จะช่วยเร่งความเร็วได้อย่างมาก โดยเฉพาะสำหรับ Model ที่ต้อง Predict ทุกๆ Bar การใช้ CUDA ผ่าน CuPy หรือ PyTorch จะช่วยให้ได้ Throughput สูงกว่า CPU หลายร้อยเท่า สำหรับการ Inference Model ขนาดเล็กบน GPU อาจเร็วกว่า CPU ถึง 1,000 เท่า ขึ้นอยู่กับ GPU Model
import cupy as cp # GPU-accelerated NumPy
import numpy as np
from numba import cuda
import time
@cuda.jit
def cuda_calculate_rsi_kernel(close_prices, periods, rsi_output, period_count):
"""
CUDA Kernel สำหรับ RSI Calculation
"""
idx = cuda.blockIdx.x * cuda.blockDim.x + cuda.threadIdx.x
if idx < len(close_prices):
for p_idx in range(period_count):
period = periods[p_idx]
if idx < period:
rsi_output[p_idx, idx] = 50.0 # Not enough data
else:
# Calculate gains and losses
gain_sum = 0.0
loss_sum = 0.0
for i in range(idx - period + 1, idx + 1):
diff = close_prices[i] - close_prices[i - 1]
if diff > 0:
gain_sum += diff
else:
loss_sum -= diff
avg_gain = gain_sum / period
avg_loss = loss_sum / period
if avg_loss == 0:
rsi_output[p_idx, idx] = 100.0
else:
rs = avg_gain / avg_loss
rsi_output[p_idx, idx] = 100.0 - (100.0 / (1.0 + rs))
def gpu_calculate_rsi(close_prices, periods):
"""Calculate RSI on GPU using CuPy"""
n = len(close_prices)
p = len(periods)
# Transfer to GPU
close_gpu = cp.asarray(close_prices)
# Allocate output on GPU
rsi_gpu = cp.zeros((p, n), dtype=cp.float64)
periods_gpu = cp.asarray(periods)
# Calculate using CuPy (simpler than writing CUDA kernels)
for i, period in enumerate(periods):
# Vectorized RSI calculation
delta = cp.diff(close_gpu, prepend=close_gpu[0])
gain = cp.maximum(delta, 0)
loss = cp.maximum(-delta, 0)
avg_gain = cp.convolve(gain, cp.ones(period)/period, mode='valid')
avg_loss = cp.convolve(loss, cp.ones(period)/period, mode='valid')
# Pad to original size
avg_gain = cp.concatenate([cp.zeros(period-1), avg_gain])
avg_loss = cp.concatenate([cp.zeros(period-1), avg_loss])
rs = avg_gain / (avg_loss + 1e-10)
rsi_gpu[i] = 100 - (100 / (1 + rs))
return cp.asnumpy(rsi_gpu)
def benchmark_gpu_vs_cpu():
"""Benchmark GPU vs CPU for RSI calculation"""
close_prices = np.random.randn(num_bars).cumsum() + 100
periods = [7, 14, 21, 28]
# CPU version
start = time.time()
cpu_result = np.zeros((len(periods), num_bars))
for i, period in enumerate(periods):
delta = np.diff(close_prices, prepend=close_prices[0])
gain = np.maximum(delta, 0)
loss = np.maximum(-delta, 0)
avg_gain = np.convolve(gain, np.ones(period)/period, mode='valid')
avg_loss = np.convolve(loss, np.ones(period)/period, mode='valid')
avg_gain = np.concatenate([np.zeros(period-1), avg_gain])
avg_loss = np.concatenate([np.zeros(period-1), avg_loss])
rs = avg_gain / (avg_loss + 1e-10)
cpu_result[i] = 100 - (100 / (1 + rs))
cpu_time = time.time() - start
# GPU version
start = time.time()
gpu_result = gpu_calculate_rsi(close_prices