Note
Go to the end to download the full example code or to run this example in your browser via JupyterLite.
Characteristics Factor Model#
This tutorial shows how to build a characteristics-based cross-sectional factor
model with CharacteristicsFactorModel. The methodology
is covered in the Factor Models user guide.
A characteristics factor model builds factor exposures from point-in-time asset characteristics (e.g. industry, market capitalization, book equity, analyst estimates). It computes the exposures from the data and estimates the factor returns at each date by cross-sectional regression of asset returns on factor exposures [1] (see Cross-Sectional Regression):
where \(B\) is the factor exposure matrix, \(f_t\) the factor returns and \(\epsilon_t\) the idiosyncratic returns.
The model then assembles the asset covariance forecast from the factor covariance \(F\), estimated on the factor return series, and the idiosyncratic covariance \(D\), estimated from the idiosyncratic returns [2]:
We build a model with 24 factors (1 global market factor, 10 industry factors and 13 style factors from 29 descriptors), fit it on a synthetic panel, run its diagnostics, evaluate its covariance forecast out of sample and tune its hyperparameters.
Data#
We use a synthetic characteristics panel generated by
make_synthetic_characteristics. It includes late
listings, delistings, holidays and missing characteristics.
Asset characteristics come from commercial point-in-time datasets whose
licences prohibit redistribution, so the gallery uses synthetic data to
present the full API. For results on real data, see the
Factor Models user guide.
import numpy as np
from plotly.io import show
from skfolio.datasets import make_synthetic_characteristics
panel = make_synthetic_characteristics(
n_assets=500, n_observations=1500, random_state=0
)
The data is stored in an AssetPanel, skfolio’s
container for aligned cross-sectional asset data. The motivation for introducing
this container and its practical benefits for portfolio workflows are covered in
the Asset Data Representation user guide.
Let’s inspect the data with to_dataframe, which converts the panel to a
pandas DataFrame:
panel.to_dataframe(output_format="long").head()
Next, we summarize the panel structure with the info method:
print(panel.info())
AssetPanel Info
============================================================
Observations : 1,500 (2015-01-01 -> 2020-09-30)
Assets : 500
Fields : 22
Panel entries : 750,000 (observations x assets)
Missing : 13.8% total, 4.5% in Active Mask
Active Mask
-----------
In mask : 677,235 / 750,000 entries (90.3%)
Assets per obs : min=429, median=456, max=467
Assets in mask : 500 / 500
median duration : 1,500 observations
shortest / longest : 264 / 1,500 observations
Estimation Mask
---------------
In mask : 645,011 / 750,000 entries (86.0%)
Assets per obs : min=410, median=434, max=445
Assets in mask : 475 / 500
median duration : 1,500 observations
shortest / longest : 264 / 1,500 observations
Field Coverage
------------------------------------------------------------
dtype % missing % missing fully missing
total in Active Mask assets (active)
returns float32 9.7% 0.0% 0
adj_close float32 9.7% 0.0% 0
adj_volume float32 9.7% 0.0% 0
adj_shares_outstanding float32 9.7% 0.0% 0
market_cap float32 9.7% 0.0% 0
ebitda_ttm float32 18.7% 10.0% 0
enterprise_value float32 18.7% 10.0% 0
net_income_ttm float32 10.6% 1.0% 0
sales_ttm float32 10.6% 1.0% 0
dividends_ttm float32 10.6% 1.0% 0
net_buybacks_ttm float32 10.6% 1.0% 0
book_equity float32 10.6% 1.0% 0
operating_cash_flow_ttm float32 10.6% 1.0% 0
total_debt float32 10.6% 1.0% 0
total_assets float32 10.6% 1.0% 0
cost_of_revenue_ttm float32 18.8% 10.0% 0
capex_ttm float32 10.6% 1.0% 0
short_interest float32 10.6% 1.0% 0
eps_ntm float32 27.8% 20.0% 0
dps_ntm float32 27.8% 20.0% 0
eps_ntm_std float32 27.7% 19.9% 0
industry int32 9.7% 0.0% 0
Categorical Fields
------------------------------------------------------------
industry : 10 levels
Min number of assets per level (over time):
< 10 : 0 levels
10 - 20 : 0 levels
20 - 50 : 8 levels (Software, Banks, Energy, Commercial, ... +4 more)
> 50 : 2 levels (Real Estate, Capital Goods)
active_mask records whether each asset belongs to the universe at each
observation, distinguishing missing data for an active asset from periods outside
the universe (e.g. before listing or after delisting). estimation_mask selects
the active asset-observation pairs used to estimate cross-sectional statistics and
factor returns. The model still computes exposures and forecasts for active pairs
outside this subset. See Coverage, Estimation and Investment Universes for details.
In this panel:
90.3% of asset-observation pairs are active.
475 of the 500 assets contribute to model estimation.
Analyst-estimate fields are missing about 20% of the time and
industryhas 10 categories.
Factor Exposures#
We define each factor with a factor exposure estimator that transforms the
panel fields into exposures (see Factor Exposures in the user guide). Estimators carry a
family attribute ("market", "industry", "style") used for
neutralization, zero-sum constraints and reporting.
Global Factor#
The global factor has unit exposure for every asset. With benchmark-weighted centering and the industry zero-sum constraint set below, its factor return captures the benchmark (market) return (see Global Factor and Benchmark Portfolio):
from skfolio.factor_exposure import GlobalFactor
global_factor = GlobalFactor(family="market")
Industry Factors#
Next, we derive the industry factor exposures by one-hot encoding the industry
categorical field:
from skfolio.factor_exposure import OneHotCategoricalFactors
industry_factors = OneHotCategoricalFactors(category="industry", family="industry")
Style Factors#
Now we build the style factors with
FixedWeightedFactor from one or more
descriptors, each computed from panel fields (e.g.
BookToPrice from book_equity and
market_cap). The estimator winsorizes and z-scores each descriptor
cross-sectionally, then combines the scores with the fixed weights
(see Descriptors and
Cross-Sectional Transformers). With
transform_by_group="industry", scoring happens within each industry,
which makes style exposures orthogonal to the industry factors (see
Neutralization):
from skfolio.descriptor import (
AssetsGrowthRate,
AssetTurnover,
BookLeverage,
BookToPrice,
CapexToAssetsChangeInIntensity,
CashFlowToAssets,
CashFlowToPrice,
DebtToAssets,
DividendToPrice,
EarningsChangeToPrice,
EarningsToPrice,
EbitdaToEnterpriseValue,
EWAmihudIlliquidity,
EWMarketBeta,
EWMomentum,
EWResidualVolatility,
EWShareTurnover,
EWVolatility,
ForwardEarningsToPrice,
GrossMargin,
GrossProfitability,
IssuanceGrowthRate,
LogMarketCap,
MarketLeverage,
ReturnOnAssets,
ReturnOnEquity,
SalesGrowthRate,
SalesToPrice,
ShareholderYield,
)
from skfolio.factor_exposure import DerivedFactor, FixedWeightedFactor
week = 5
month = 21
quarter = 3 * month
half_year = 6 * month
year = 12 * month
beta_factor = FixedWeightedFactor(
descriptors=[("market_beta", EWMarketBeta(half_life=year))],
transform_by_group="industry",
)
momentum_factor = FixedWeightedFactor(
descriptors=[("momentum", EWMomentum(half_life=half_year, skip=month))],
transform_by_group="industry",
)
size_factor = FixedWeightedFactor(
descriptors=[("log_mcap", LogMarketCap())], transform_by_group="industry"
)
non_linear_size_factor = DerivedFactor(
source="size", func=lambda x: x**3, transform_by_group="industry"
)
value_factor = FixedWeightedFactor(
descriptors=[
("book_to_price", BookToPrice()),
("sales_to_price", SalesToPrice()),
("cash_flow_to_price", CashFlowToPrice()),
],
weights=[0.8, 0.1, 0.1],
transform_by_group="industry",
)
earnings_yield_factor = FixedWeightedFactor(
descriptors=[
("fwd_earnings_to_price", ForwardEarningsToPrice()),
("earnings_to_price", EarningsToPrice()),
("enterprise_multiple", EbitdaToEnterpriseValue()),
],
transform_by_group="industry",
)
growth_factor = FixedWeightedFactor(
descriptors=[
("earnings_change_to_price", EarningsChangeToPrice(lag=year)),
("sales_growth", SalesGrowthRate(lag=year)),
],
transform_by_group="industry",
)
profitability_factor = FixedWeightedFactor(
descriptors=[
("asset_turnover", AssetTurnover()),
("gross_profitability", GrossProfitability()),
("gross_margin", GrossMargin()),
("return_on_assets", ReturnOnAssets()),
("return_on_equity", ReturnOnEquity()),
("cash_flow_to_assets", CashFlowToAssets()),
],
transform_by_group="industry",
)
investment_factor = FixedWeightedFactor(
descriptors=[
("asset_growth", AssetsGrowthRate(lag=year)),
("issuance_growth", IssuanceGrowthRate(lag=year)),
("capex_growth", CapexToAssetsChangeInIntensity(lag=year)),
],
transform_by_group="industry",
)
dividend_yield_factor = FixedWeightedFactor(
descriptors=[
("dividend_to_price", DividendToPrice()),
("shareholder_yield", ShareholderYield()),
],
weights=[0.7, 0.3],
transform_by_group="industry",
)
leverage_factor = FixedWeightedFactor(
descriptors=[
("market_leverage", MarketLeverage()),
("debt_to_assets", DebtToAssets()),
("book_leverage", BookLeverage()),
],
transform_by_group="industry",
)
liquidity_factor = FixedWeightedFactor(
descriptors=[
("share_turnover", EWShareTurnover(half_life=quarter)),
("amihud_illiquidity", EWAmihudIlliquidity()),
],
transform_by_group="industry",
)
volatility_factor = FixedWeightedFactor(
descriptors=[
("vol", EWVolatility(half_life=quarter)),
(
"residual_vol",
EWResidualVolatility(half_life=quarter, beta_half_life=quarter),
),
],
transform_by_group="industry",
)
The non-linear size factor is a
DerivedFactor computed as the cube of the
size exposure, capturing return patterns that differ between mid-caps and
the extremes of the size spectrum. The momentum descriptor
EWMomentum skips the most recent month
(skip=month) to separate medium-term momentum from short-term reversal.
Model Definition#
We now assemble the factors into the
CharacteristicsFactorModel:
from skfolio.moments import EWMu, RegimeAdjustedEWCovariance
from skfolio.prior import CharacteristicsFactorModel, EmpiricalPrior
model = CharacteristicsFactorModel(
factors=[
("market", global_factor),
("industry", industry_factors),
("beta", beta_factor),
("momentum", momentum_factor),
("size", size_factor),
("non_linear_size", non_linear_size_factor),
("value", value_factor),
("earnings_yield", earnings_yield_factor),
("growth", growth_factor),
("profitability", profitability_factor),
("investment", investment_factor),
("dividend_yield", dividend_yield_factor),
("leverage", leverage_factor),
("liquidity", liquidity_factor),
("volatility", volatility_factor),
],
neutralize_against={
"non_linear_size": ["size"],
"volatility": ["beta"],
},
constrained_families=[("industry", None)],
exposure_lag=1,
inv_idio_variance_weight_shrinkage=0.5,
factor_prior_estimator=EmpiricalPrior(
covariance_estimator=RegimeAdjustedEWCovariance(
half_life=half_year, corr_half_life=year, regime_half_life=month
),
mu_estimator=EWMu(half_life=year),
),
n_jobs=-1,
)
Let’s walk through the main parameters:
neutralize_againstorthogonalizes the non-linear size exposure against size and the volatility exposure against beta, keeping only the component uncorrelated with the target factor (see Neutralization).constrained_familiesapplies a benchmark-weighted zero-sum constraint on the industry factor returns. The one-hot industry exposures sum to the market exposure, so without a constraint the regression design is rank-deficient (see Zero-Sum Constraints). With it, the market factor captures the benchmark return and industry factors capture relative effects around it.exposure_lag=1regresses the returns over \((t-1, t]\) on exposures measured at \(t-1\), avoiding look-ahead bias (see Time Alignment and Look-Ahead Bias).inv_idio_variance_weight_shrinkage=0.5blends square-root market-cap regression weights with inverse-idiosyncratic-variance weights in a two-pass feasible GLS (see Regression Weights).The factor prior estimates expected factor returns with an exponentially weighted mean and the factor covariance with a regime-adjusted exponentially weighted estimator (see Factor Return Distribution).
Fitting#
Full-Panel Fitting#
We fit the model on the panel. The optional X argument defines the
investment universe used by downstream portfolio optimization. When omitted,
it defaults to the panel’s coverage universe. See Coverage, Estimation
and Investment Universes.
model.fit(characteristics=panel)
Online Learning#
The model also supports online learning (see Online
Learning). partial_fit appends new
observations without refitting the history, and the result is identical to
a full-panel fit on the concatenated data. We initialize the model with two
years plus one month of data to cover the cumulative descriptor and estimator
warmups (see Warmup Periods). Later chunks can
be as small as one observation. The equivalent incremental fit is:
warmup = 2 * year + month
model.fit(characteristics=panel[:warmup])
for i in range(warmup, len(panel), month):
model.partial_fit(characteristics=panel[i : i + month])
Model Outputs#
Now let’s look at what the fitted model produces. return_distribution_
holds the ReturnDistribution consumed by all
skfolio optimizations: expected returns mu, asset covariance covariance
and asset return scenarios returns. Assets that are not investable at the last
observation (e.g. delisted or still in warmup) carry NaN moments. Compatible
optimizers solve on the investable subset and assign zero weight to the
rest:
distribution = model.return_distribution_
print(f"mu shape: {distribution.mu.shape}")
print(f"covariance shape: {distribution.covariance.shape}")
print(f"scenarios shape: {distribution.returns.shape}")
print(f"investable assets: {np.isfinite(distribution.mu).sum()}")
mu shape: (500,)
covariance shape: (500, 500)
scenarios shape: (1247, 500)
investable assets: 418
The scenario history contains 1,247 observations because the one-year descriptor warmup and one-observation exposure lag consume the first 253 observations of the 1,500-observation panel.
Next, we retrieve factor_model_, the FactorModel
container holding the full decomposition: factor exposures, loading matrix,
factor returns, expected factor returns and covariance, idiosyncratic
returns, variances and covariance, regression and benchmark weights, plus
the diagnostics used in the sections below. Its summary method reports
per-factor statistics:
factor_model = model.factor_model_
factor_model.summary(families=["market", "style"])
These are pure-factor returns. Each factor return is the cross-sectional regression coefficient for one unit of exposure, the return of a factor-mimicking portfolio with unit exposure to that factor and zero exposure to all others.
Exposure Diagnostics#
Now we check the conditioning of the exposure matrix (see Exposure Diagnostics). Strong collinearity between exposures inflates the variance of the estimated factor returns and makes attribution unstable. We start with the time-average pairwise exposure correlations:
fig = factor_model.plot_exposure_correlation()
fig.update_layout(height=700)
The market exposure is constant (every asset has unit exposure), so its
correlation row is uninformative and displays as zero. Market neutrality is
instead guaranteed by benchmark-weighted centering: every exposure has a
zero benchmark-weighted mean, so no factor carries net market exposure.
The style-industry correlations are zero, the result of
within-industry scoring
(transform_by_group="industry"), and so are the volatility-beta and
non-linear size-size correlations, the result of neutralize_against.
Industry factors are slightly negatively correlated rather than
uncorrelated. Each asset belongs to exactly one industry, so membership in
one industry rules out membership in all others. Zero correlation would
mean industry memberships are independent. The remaining style
correlations reflect the correlated traits of the synthetic universe (e.g.
value with leverage and earnings yield), with no redundant factors,
consistent with the low mean variance inflation factors (mean_vif)
reported in the summary table above.
Note
For two one-hot exposures \(x\) and \(y\) with benchmark weights \(p_x\) and \(p_y\), the product \(xy\) is always zero, so the covariance \(\mathbb{E}[xy] - \mathbb{E}[x]\mathbb{E}[y] = -p_x p_y\) is negative, giving a correlation of \(-\sqrt{p_x p_y / ((1-p_x)(1-p_y))}\), about -0.11 for ten industries of similar weight.
Regression Diagnostics#
Let’s inspect the fit of the cross-sectional regressions with
cs_regression_scores, which reports per-observation fit statistics:
factor_model.cs_regression_scores.mean()
r2 0.608784
adjusted_r2 0.585589
aic -4013.061835
bic -3920.574305
dtype: float64
factor_model.plot_cs_regression_scores(score="adjusted_r2", window=20)
The mean \(R^2\) is around 60%, as expected for this synthetic dataset, where factor effects are deliberately easy to identify. On real daily US equity data, the mean \(R^2\) typically falls between 25% and 40% (see Regression Diagnostics in the user guide). The synthetic universe is intentionally smaller and cleaner than real-world datasets.
Next, plot_cs_regression_t_stat_exceedance_rate shows how often each factor’s
t-statistic exceeds 2 in absolute value. A factor whose true coefficient is
zero would exceed it about 5% of the time, so rates persistently above this
reference identify factors that are repeatedly significant in the
cross-section:
factor_model.plot_cs_regression_t_stat_exceedance_rate(families=["market", "style"])
Factor Returns#
Let’s plot the estimated factor returns accumulated through time:
fig = factor_model.plot_factor_cumulative_returns(families=["market", "style"])
show(fig)
The market factor tracks the cumulative benchmark return, as expected from the benchmark-weighted centering and the industry zero-sum constraint. Among the style factors, momentum shows the largest cumulative return and the highest annualized Sharpe ratio in the summary table above. The synthetic generator assigns momentum a persistent premium, and the model recovers it from the cross-section of returns.
Idiosyncratic Risk Calibration#
Now let’s test the idiosyncratic volatility forecasts. The diagnostics use the standardized idiosyncratic returns \(z_{it} = \epsilon_{it} / \hat\sigma_{it}\). Under correct calibration, \(z\) has cross-sectional standard deviation 1.0:
factor_model.idio_calibration_summary()
mean_cs_std 1.001165
median_cs_std 0.998134
mean_cs_excess_kurtosis 0.352186
mean_cs_skewness 0.003433
mean_tail_rate_3sigma 0.004793
Name: idio_calibration, dtype: float64
factor_model.plot_idio_calibration(window=20)
The series oscillates around 1.0, so the idiosyncratic risk forecasts are
correctly scaled. The mean_tail_rate_3sigma of about 0.5%, above the
0.27% Gaussian reference, reflects the fat tails of the idiosyncratic
returns, a feature of real equity data reproduced by the generator through
Student-t shocks.
Let’s now separate ranking power from calibration with plot_idio_vol_ic,
which displays the Spearman correlation between idiosyncratic volatility
forecasts and next-period absolute idiosyncratic returns. High values mean
the model tends to forecast higher volatility for assets that subsequently
experience larger idiosyncratic moves. This tests the ordering of the
forecasts, while the calibration series above tests their scale. Together,
they show that the idiosyncratic risk forecasts are both well ordered and
well scaled:
factor_model.plot_idio_vol_ic()
Information Coefficient#
Finally, we measure the return-predictive power of the exposures with
exposure_ic_summary, the cross-sectional rank correlation between factor
exposures at \(t\) and forward asset returns:
factor_model.exposure_ic_summary(families=["market", "style"])
factor_model.plot_cumulative_exposure_ic(families=["market", "style"])
Mean ICs are small, IC IRs remain low and hit rates remain around 50% across factors, so no factor stands out as a strong, consistent predictor of next-period returns.
The IC quantifies return-predictive power [3]. In a risk model, factors are designed to forecast covariance, not expected returns. A factor can therefore be an excellent risk factor with an IC near zero, and a low IC is not a reason to discard it. The IC is mainly useful for evaluating alpha signals and factor premia, not for deciding whether a factor should remain in a risk model.
Covariance Forecast Evaluation#
We now evaluate the full covariance forecast out of sample with
online_covariance_forecast_evaluation (see
Covariance Forecast Evaluation). It walks forward through
the data, compares each covariance forecast with the subsequently realized
returns and computes calibration diagnostics.
X are the asset returns in the standard skfolio format and
warmup_size reserves two years plus one month for the estimator warmups
(see Warmup Periods):
from skfolio.model_selection import online_covariance_forecast_evaluation
X = panel.to_dataframe(fields="returns")
evaluation = online_covariance_forecast_evaluation(
model,
X,
params={"characteristics": panel},
warmup_size= 2 * year + month,
test_size=week,
)
evaluation.summary()
The Mahalanobis ratio tests the full covariance structure and the diagonal
ratio the individual asset variances, both with a 1.0 target. Values above
1.0 indicate underestimated risk and below 1.0 overestimated risk. The
portfolio standardized returns test calibration along a portfolio direction,
with std as the bias statistic, and the portfolio QLIKE scores the
portfolio variance forecasts, lower being better.
evaluation.plot_calibration(diagnostics=["bias"])
The bias statistic oscillates around its 1.0 target, and the summary table shows the Mahalanobis and diagonal calibration ratios close to 1.0 as well, so the risk forecasts are well calibrated across the sample.
Hyperparameter Tuning#
Model parameters follow the scikit-learn naming convention. Here we jointly tune the regression weights, factor covariance dynamics and idiosyncratic variance dynamics.
Because these parameters are continuous, we use
OnlineRandomizedSearch. It evaluates each
candidate in a single walk-forward pass using partial_fit. The search uses
the first four years, preserving the remaining observations as an untouched
holdout.
We score the risk model directly with a portfolio-variance QLIKE loss. In
make_scorer, response_method=None indicates a non-predictor estimator and
greater_is_better=False a loss to minimize. The complete search is shown
without execution because fitting 100 candidates is computationally
intensive (see Hyperparameter Tuning):
from scipy.stats import loguniform, uniform
from skfolio.metrics import make_scorer, portfolio_variance_qlike_loss
from skfolio.model_selection import OnlineRandomizedSearch
qlike_scorer = make_scorer(
portfolio_variance_qlike_loss,
greater_is_better=False,
response_method=None,
)
search_size = 4 * year
search = OnlineRandomizedSearch(
estimator=model,
param_distributions={
"inv_idio_variance_weight_shrinkage": uniform(0.0, 1.0),
"factor_prior_estimator__covariance_estimator__half_life": (
loguniform(quarter, half_year)
),
"factor_prior_estimator__covariance_estimator__corr_half_life": (
loguniform(half_year, year)
),
"factor_prior_estimator__covariance_estimator__regime_half_life": (
loguniform(week, quarter)
),
"idio_variance_estimator__half_life": loguniform(
month, half_year
),
"idio_variance_estimator__regime_half_life": loguniform(
week, quarter
),
},
n_iter=100,
scoring=qlike_scorer,
warmup_size=2 * year + month,
test_size=week,
random_state=0,
n_jobs=-1,
)
search.fit(
X.iloc[:search_size],
characteristics=panel[:search_size],
)
print(search.best_params_)
print(search.best_score_)
best_score_ summarizes the walk-forward validation windows within the
four-year search period. cv_results_ contains every sampled parameter set,
score, rank and fit time. For a small set of discrete candidates,
OnlineGridSearch is also available and
evaluates every parameter combination.
Conclusion#
We built a 24-factor characteristics factor model, verified its exposures
and regressions, evaluated its risk forecasts out of sample, showed how to
tune them and updated the model online. The fitted model is a prior estimator,
so we can pass it to any skfolio optimizer through prior_estimator. It
supplies expected returns, covariance, scenarios and factor exposures that
can be constrained.
The next tutorial builds a factor-constrained portfolio on top of this model and performs factor-level risk and performance attribution.
See also
The Factor Models user guide covers the methodology in depth, including portfolio construction, attribution and alpha integration with this model.
References#
Total running time of the script: (0 minutes 44.412 seconds)
