AAPL classification from SPY regimes

This example uses homegrown q.indicator functions to turn SPY trend into predefined short, hold, and long targets. A scikit-learn model then learns to predict the next session’s SPY regime from AAPL market features.

The target is an interpretable indicator rule rather than a claim about future returns. Treat the workflow as an evaluation template, not as a trading strategy.

import numpy as np
import pandas as pd
import plotly.io as pio
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import log_loss
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

import qrt as q

pio.renderers.default = "notebook_connected"

1. Get data

Load the bundled AAPL and SPY histories and align their trading sessions. The model features describe AAPL using returns, trend, volatility, range, volume, and AAPL’s relative strength against SPY.

aapl = q.data.datasets.load("aapl")
spy = q.data.datasets.load("spy")
common_sessions = aapl.index.intersection(spy.index)
aapl = aapl.loc[common_sessions]
spy = spy.loc[common_sessions]

aapl_returns = aapl["close"].pct_change()
relative_strength = q.indicator.relative_strength(
    aapl["close"], spy["close"], lookback=20
)
relative_strength_average = q.indicator.rsma(relative_strength, window=20)
relative_strength_phase = q.indicator.rs_phase(
    relative_strength, relative_strength_average
)
features = pd.DataFrame(
    {
        "return_1d": aapl_returns,
        "momentum_5d": aapl["close"].pct_change(5),
        "momentum_20d": aapl["close"].pct_change(20),
        "close_vs_sma_20d": aapl["close"].div(
            q.indicator.sma(aapl["close"], 20)
        ).sub(1),
        "ema_20d_vs_sma_100d": q.indicator.ema(aapl["close"], 20).div(
            q.indicator.sma(aapl["close"], 100)
        ).sub(1),
        "relative_strength_20d": relative_strength,
        "rsma_20d": relative_strength_average,
        "rs_phase": relative_strength_phase["rs_phase"].astype(int),
        "rs_phase_days": relative_strength_phase["rs_phase_days"],
        "rs_days_in_correction_60d": q.indicator.rs_days(
            relative_strength, spy["close"], window=60
        ),
        "rs_new_high_before_price_60d": q.indicator.rsnhbp(
            aapl["close"], relative_strength, window=60
        ).astype(int),
        "volatility_20d": aapl_returns.rolling(20).std(),
        "range_1d": (aapl["high"] - aapl["low"]) / aapl["close"],
        "volume_change_5d": aapl["volume"].pct_change(5),
    }
).replace([np.inf, -np.inf], np.nan).dropna()

features.tail()

	return_1d	momentum_5d	momentum_20d	close_vs_sma_20d	ema_20d_vs_sma_100d	relative_strength_20d	rsma_20d	rs_phase	rs_phase_days	rs_days_in_correction_60d	rs_new_high_before_price_60d	volatility_20d	range_1d	volume_change_5d
datetime
2026-07-17	0.001440	0.058417	0.127690	0.092374	0.104660	0.121962	0.037196	1	11	0.0	0	0.022858	0.017948	0.857686
2026-07-20	-0.021424	0.029246	0.095903	0.063995	0.107562	0.102130	0.043380	1	12	0.0	0	0.023679	0.030711	0.236031
2026-07-21	0.003521	0.040907	0.103464	0.062423	0.110397	0.098239	0.048605	1	13	0.0	0	0.023603	0.022518	0.137659
2026-07-22	-0.005645	-0.004916	0.107340	0.051045	0.112154	0.088487	0.052403	1	14	0.0	0	0.023505	0.017368	-0.364215
2026-07-23	-0.012980	-0.034808	0.097516	0.032643	0.111944	0.090779	0.056058	1	15	0.0	0	0.023775	0.012280	-0.351431

2. Split chronologically

Use the oldest 60% of AAPL observations for training, the next 20% for validation, and the newest 20% for testing. One session is purged before each boundary because every feature row will receive the following session’s regime label.

train_boundary = int(len(features) * 0.60)
validation_boundary = int(len(features) * 0.80)

X = {
    "train": features.iloc[: train_boundary - 1],
    "validation": features.iloc[train_boundary : validation_boundary - 1],
    "test": features.iloc[validation_boundary:-1],
}

pd.DataFrame(
    {
        name: {
            "start": frame.index.min().date(),
            "end": frame.index.max().date(),
            "observations": len(frame),
        }
        for name, frame in X.items()
    }
).T

	start	end	observations
train	2000-05-24	2016-02-01	3946
validation	2016-02-03	2021-04-23	1315
test	2021-04-27	2026-07-22	1315

3. Label with SPY indicators

Define the target from SPY’s q.indicator.ema(20) relative to its q.indicator.sma(100). A spread above 2% is long, below -2% is short, and the neutral band is hold. These fixed rules are applied unchanged to train, validation, and test.

classes = ["short", "hold", "long"]
regime_band = 0.02
spy_trend_spread = q.indicator.ema(spy["close"], 20).div(
    q.indicator.sma(spy["close"], 100)
).sub(1)
spy_regime = pd.Series(
    np.select(
        [spy_trend_spread < -regime_band, spy_trend_spread > regime_band],
        ["short", "long"],
        default="hold",
    ),
    index=spy.index,
    name="position",
)
next_session_regime = spy_regime.shift(-1)
y = {
    name: next_session_regime.loc[frame.index]
    for name, frame in X.items()
}

pd.DataFrame(
    {name: labels.value_counts(normalize=True) for name, labels in y.items()}
).T.loc[:, classes].style.format("{:.1%}")

position	short	hold	long
train	20.5%	35.6%	44.0%
validation	10.6%	22.3%	67.1%
test	15.1%	21.6%	63.3%

4. Train with scikit-learn

Build a scikit-learn pipeline that standardizes the AAPL features and fits multinomial logistic regression. Keeping preprocessing inside the pipeline ensures its statistics are learned from training data only.

def make_model(c_value):
    return make_pipeline(
        StandardScaler(),
        LogisticRegression(C=c_value, max_iter=2_000),
    )


model = make_model(c_value=1.0)
model.fit(X["train"], y["train"])
model

Pipeline(steps=[('standardscaler', StandardScaler()),
                ('logisticregression', LogisticRegression(max_iter=2000))])

In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.

Pipeline

?Documentation for PipelineiFitted

Parameters

	steps steps: list of tuples List of (name of step, estimator) tuples that are to be chained in sequential order. To be compatible with the scikit-learn API, all steps must define `fit`. All non-last steps must also define `transform`. See :ref:`Combining Estimators <combining_estimators>` for more details.	[('standardscaler', ...), ('logisticregression', ...)]
	transform_input transform_input: list of str, default=None The names of the :term:`metadata` parameters that should be transformed by the pipeline before passing it to the step consuming it. This enables transforming some input arguments to ``fit`` (other than ``X``) to be transformed by the steps of the pipeline up to the step which requires them. Requirement is defined via :ref:`metadata routing <metadata_routing>`. For instance, this can be used to pass a validation set through the pipeline. You can only set this if metadata routing is enabled, which you can enable using ``sklearn.set_config(enable_metadata_routing=True)``. .. versionadded:: 1.6	None
	memory memory: str or object with the joblib.Memory interface, default=None Used to cache the fitted transformers of the pipeline. The last step will never be cached, even if it is a transformer. By default, no caching is performed. If a string is given, it is the path to the caching directory. Enabling caching triggers a clone of the transformers before fitting. Therefore, the transformer instance given to the pipeline cannot be inspected directly. Use the attribute ``named_steps`` or ``steps`` to inspect estimators within the pipeline. Caching the transformers is advantageous when fitting is time consuming. See :ref:`sphx_glr_auto_examples_neighbors_plot_caching_nearest_neighbors.py` for an example on how to enable caching.	None
	verbose verbose: bool, default=False If True, the time elapsed while fitting each step will be printed as it is completed.	False

Fitted attributes

Name	Type	Value
classes_ classes_: ndarray of shape (n_classes,) The classes labels. Only exist if the last step of the pipeline is a classifier.	ndarray[object](3,)	['hold','long','short']
feature_names_in_ feature_names_in_: ndarray of shape (`n_features_in_`,) Names of features seen during :term:`fit`. Only defined if the underlying estimator exposes such an attribute when fit. .. versionadded:: 1.0	ndarray[object](14,)	['return_1d','momentum_5d','momentum_20d',...,'volatility_20d','range_1d', 'volume_change_5d']
n_features_in_ n_features_in_: int Number of features seen during :term:`fit`. Only defined if the underlying first estimator in `steps` exposes such an attribute when fit. .. versionadded:: 0.24	int	14

StandardScaler

?Documentation for StandardScaler

Parameters

	copy copy: bool, default=True If False, try to avoid a copy and do inplace scaling instead. This is not guaranteed to always work inplace; e.g. if the data is not a NumPy array or scipy.sparse CSR matrix, a copy may still be returned.	True
	with_mean with_mean: bool, default=True If True, center the data before scaling. This does not work (and will raise an exception) when attempted on sparse matrices, because centering them entails building a dense matrix which in common use cases is likely to be too large to fit in memory.	True
	with_std with_std: bool, default=True If True, scale the data to unit variance (or equivalently, unit standard deviation).	True

Fitted attributes

Name	Type	Value
feature_names_in_ feature_names_in_: ndarray of shape (`n_features_in_`,) Names of features seen during :term:`fit`. Defined only when `X` has feature names that are all strings. .. versionadded:: 1.0	ndarray[object](14,)	['return_1d','momentum_5d','momentum_20d',...,'volatility_20d','range_1d', 'volume_change_5d']
mean_ mean_: ndarray of shape (n_features,) or None The mean value for each feature in the training set. Equal to ``None`` when ``with_mean=False`` and ``with_std=False``.	ndarray[float64](14,)	[0. ,0.01,0.02,...,0.02,0.03,0.14]
n_features_in_ n_features_in_: int Number of features seen during :term:`fit`. .. versionadded:: 0.24	int	14
n_samples_seen_ n_samples_seen_: int or ndarray of shape (n_features,) The number of samples processed by the estimator for each feature. If there are no missing samples, the ``n_samples_seen`` will be an integer, otherwise it will be an array of dtype int. If `sample_weights` are used it will be a float (if no missing data) or an array of dtype float that sums the weights seen so far. Will be reset on new calls to fit, but increments across ``partial_fit`` calls.	float64	3946
scale_ scale_: ndarray of shape (n_features,) or None Per feature relative scaling of the data to achieve zero mean and unit variance. Generally this is calculated using `np.sqrt(var_)`. If a variance is zero, we can't achieve unit variance, and the data is left as-is, giving a scaling factor of 1. `scale_` is equal to `None` when `with_std=False`. .. versionadded:: 0.17 scale_	ndarray[float64](14,)	[0.03,0.06,0.12,...,0.01,0.02,0.72]
var_ var_: ndarray of shape (n_features,) or None The variance for each feature in the training set. Used to compute `scale_`. Equal to ``None`` when ``with_mean=False`` and ``with_std=False``.	ndarray[float64](14,)	[0. ,0. ,0.01,...,0. ,0. ,0.52]

14 features

return_1d
momentum_5d
momentum_20d
close_vs_sma_20d
ema_20d_vs_sma_100d
relative_strength_20d
rsma_20d
rs_phase
rs_phase_days
rs_days_in_correction_60d
rs_new_high_before_price_60d
volatility_20d
range_1d
volume_change_5d

LogisticRegression

?Documentation for LogisticRegression

Parameters

	max_iter max_iter: int, default=100 Maximum number of iterations taken for the solvers to converge.	2000
	penalty penalty: {'l1', 'l2', 'elasticnet', None}, default='l2' Specify the norm of the penalty: - `None`: no penalty is added; - `'l2'`: add an L2 penalty term and it is the default choice; - `'l1'`: add an L1 penalty term; - `'elasticnet'`: both L1 and L2 penalty terms are added. .. warning:: Some penalties may not work with some solvers. See the parameter `solver` below, to know the compatibility between the penalty and solver. .. versionadded:: 0.19 l1 penalty with SAGA solver (allowing 'multinomial' + L1) .. deprecated:: 1.8 `penalty` was deprecated in version 1.8 and will be removed in 1.10. Use `l1_ratio` and `C` instead. `l1_ratio=0` for `penalty='l2'`, `l1_ratio=1` for `penalty='l1'`, `l1_ratio` set to any float between 0 and 1 for `penalty='elasticnet'`, and `C=np.inf` for `penalty=None`.	'deprecated'
	C C: float, default=1.0 Inverse of regularization strength; must be a positive float. Like in support vector machines, smaller values specify stronger regularization. `C=np.inf` results in unpenalized logistic regression. For a visual example on the effect of tuning the `C` parameter with an L1 penalty, see: :ref:`sphx_glr_auto_examples_linear_model_plot_logistic_path.py`.	1.0
	l1_ratio l1_ratio: float, default=0.0 The Elastic-Net mixing parameter, with `0 <= l1_ratio <= 1`. Setting `l1_ratio=1` gives a pure L1-penalty, setting `l1_ratio=0` a pure L2-penalty. Any value between 0 and 1 gives an Elastic-Net penalty of the form `l1_ratio * L1 + (1 - l1_ratio) * L2`. .. warning:: Certain values of `l1_ratio`, i.e. some penalties, may not work with some solvers. See the parameter `solver` below, to know the compatibility between the penalty and solver. .. versionchanged:: 1.8 Default value changed from None to 0.0. .. deprecated:: 1.8 `None` is deprecated and will be removed in version 1.10. Always use `l1_ratio` to specify the penalty type.	0.0
	dual dual: bool, default=False Dual (constrained) or primal (regularized, see also :ref:`this equation <regularized-logistic-loss>`) formulation. Dual formulation is only implemented for l2 penalty with liblinear solver. Prefer `dual=False` when n_samples > n_features.	False
	tol tol: float, default=1e-4 Tolerance for stopping criteria.	0.0001
	fit_intercept fit_intercept: bool, default=True Specifies if a constant (a.k.a. bias or intercept) should be added to the decision function.	True
	intercept_scaling intercept_scaling: float, default=1 Useful only when the solver `liblinear` is used and `self.fit_intercept` is set to `True`. In this case, `x` becomes `[x, self.intercept_scaling]`, i.e. a "synthetic" feature with constant value equal to `intercept_scaling` is appended to the instance vector. The intercept becomes ``intercept_scaling * synthetic_feature_weight``. .. note:: The synthetic feature weight is subject to L1 or L2 regularization as all other features. To lessen the effect of regularization on synthetic feature weight (and therefore on the intercept) `intercept_scaling` has to be increased.	1
	class_weight class_weight: dict or 'balanced', default=None Weights associated with classes in the form ``{class_label: weight}``. If not given, all classes are supposed to have weight one. The "balanced" mode uses the values of y to automatically adjust weights inversely proportional to class frequencies in the input data as ``n_samples / (n_classes * np.bincount(y))``. Note that these weights will be multiplied with sample_weight (passed through the fit method) if sample_weight is specified. .. versionadded:: 0.17 class_weight='balanced'	None
	random_state random_state: int, RandomState instance, default=None Used when ``solver`` == 'sag', 'saga' or 'liblinear' to shuffle the data. See :term:`Glossary <random_state>` for details.	None
	solver solver: {'lbfgs', 'liblinear', 'newton-cg', 'newton-cholesky', 'sag', 'saga'}, default='lbfgs' Algorithm to use in the optimization problem. Default is 'lbfgs'. To choose a solver, you might want to consider the following aspects: - 'lbfgs' is a good default solver because it works reasonably well for a wide class of problems. - For :term:`multiclass` problems (`n_classes >= 3`), all solvers except 'liblinear' minimize the full multinomial loss, 'liblinear' will raise an error. - 'newton-cholesky' is a good choice for `n_samples` >> `n_features * n_classes`, especially with one-hot encoded categorical features with rare categories. Be aware that the memory usage of this solver has a quadratic dependency on `n_features * n_classes` because it explicitly computes the full Hessian matrix. - For small datasets, 'liblinear' is a good choice, whereas 'sag' and 'saga' are faster for large ones; - 'liblinear' can only handle binary classification by default. To apply a one-versus-rest scheme for the multiclass setting one can wrap it with the :class:`~sklearn.multiclass.OneVsRestClassifier`. .. warning:: The choice of the algorithm depends on the penalty chosen (`l1_ratio=0` for L2-penalty, `l1_ratio=1` for L1-penalty and `0 < l1_ratio < 1` for Elastic-Net) and on (multinomial) multiclass support: ================= ======================== ====================== solver l1_ratio multinomial multiclass ================= ======================== ====================== 'lbfgs' l1_ratio=0 yes 'liblinear' l1_ratio=1 or l1_ratio=0 no 'newton-cg' l1_ratio=0 yes 'newton-cholesky' l1_ratio=0 yes 'sag' l1_ratio=0 yes 'saga' 0<=l1_ratio<=1 yes ================= ======================== ====================== .. note:: 'sag' and 'saga' fast convergence is only guaranteed on features with approximately the same scale. You can preprocess the data with a scaler from :mod:`sklearn.preprocessing`. .. seealso:: Refer to the :ref:`User Guide <Logistic_regression>` for more information regarding :class:`LogisticRegression` and more specifically the :ref:`Table <logistic_regression_solvers>` summarizing solver/penalty supports. .. versionadded:: 0.17 Stochastic Average Gradient (SAG) descent solver. Multinomial support in version 0.18. .. versionadded:: 0.19 SAGA solver. .. versionchanged:: 0.22 The default solver changed from 'liblinear' to 'lbfgs' in 0.22. .. versionadded:: 1.2 newton-cholesky solver. Multinomial support in version 1.6.	'lbfgs'
	verbose verbose: int, default=0 For the liblinear and lbfgs solvers set verbose to any positive number for verbosity.	0
	warm_start warm_start: bool, default=False When set to True, reuse the solution of the previous call to fit as initialization, otherwise, just erase the previous solution. Useless for liblinear solver. See :term:`the Glossary <warm_start>`. .. versionadded:: 0.17 warm_start to support lbfgs, newton-cg, sag, saga solvers.	False
	n_jobs n_jobs: int, default=None Does not have any effect. .. deprecated:: 1.8 `n_jobs` is deprecated in version 1.8 and will be removed in 1.10.	None

Fitted attributes

Name	Type	Value
classes_ classes_: ndarray of shape (n_classes, ) A list of class labels known to the classifier.	ndarray[object](3,)	['hold','long','short']
coef_ coef_: ndarray or CSR matrix of shape (1, n_features) or (n_classes, n_features) Coefficients of the features in the decision function. `coef_` is of shape (1, n_features) when the given problem is binary. By default, it will be created as a dense array, but can be turned to sparse (CSR format) through :meth:`sparsify` (which can be beneficial under L1 regularization when many coefficients are zero), and back to dense through :meth:`densify`.	ndarray[float64](3, 14)	[[ 0. ,-0.03,-0.29,..., 0.31,-0.08, 0. ], [-0.01,-0.07, 3.18,...,-0.66, 0.07, 0.08], [ 0.01, 0.1 ,-2.89,..., 0.35, 0. ,-0.08]]
intercept_ intercept_: ndarray of shape (1,) or (n_classes,) Intercept (a.k.a. bias) added to the decision function. If `fit_intercept` is set to False, the intercept is set to zero. `intercept_` is of shape (1,) when the given problem is binary.	ndarray[float64](3,)	[ 1.14, 0.2 ,-1.34]
n_features_in_ n_features_in_: int Number of features seen during :term:`fit`. .. versionadded:: 0.24	int	14
n_iter_ n_iter_: ndarray of shape (1, ) Actual number of iterations for all classes. .. versionchanged:: 0.20 In SciPy <= 1.0.0 the number of lbfgs iterations may exceed ``max_iter``. ``n_iter_`` will now report at most ``max_iter``.	ndarray[int32](1,)	[40]

5. Validate

Fit a small inverse-regularization grid and select the value with the lowest scikit-learn log loss on validation probabilities. The test partition remains untouched during this choice.

models = {}
validation_rows = []
for c_value in [0.01, 0.1, 1.0, 10.0, 100.0]:
    candidate = make_model(c_value)
    candidate.fit(X["train"], y["train"])
    validation_score = pd.DataFrame(
        candidate.predict_proba(X["validation"]),
        index=X["validation"].index,
        columns=candidate.classes_,
    )
    models[c_value] = candidate
    validation_rows.append(
        {
            "C": c_value,
            "log_loss": log_loss(
                y["validation"],
                validation_score,
                labels=validation_score.columns,
            ),
        }
    )

validation_results = pd.DataFrame(validation_rows)
best_c = validation_results.loc[validation_results["log_loss"].idxmin(), "C"]
model = models[best_c]
validation_results.style.format({"C": "{:.3g}", "log_loss": "{:.3f}"})

	C	log_loss
0	0.01	0.588
1	0.1	0.518
2	1	0.508
3	10	0.508
4	100	0.508

6. Test

Call predict_proba once on the held-out AAPL test period, then pass the resulting class probabilities to q.stats. The tidy one-vs-rest curves remain available independently of Plotly.

y_true = y["test"]
y_score = pd.DataFrame(
    model.predict_proba(X["test"]),
    index=X["test"].index,
    columns=model.classes_,
).loc[:, classes]

roc_curves = q.stats.multiclass_roc_curve(y_true, y_score)
pr_curves = q.stats.multiclass_precision_recall_curve(y_true, y_score)

roc_summary = roc_curves.assign(
    label=lambda frame: frame["class"].fillna(frame["curve"])
).groupby("label", sort=False)["auc"].first()
pr_summary = pr_curves.assign(
    label=lambda frame: frame["class"].fillna(frame["curve"])
).groupby("label", sort=False)["average_precision"].first()

pd.concat(
    [roc_summary.rename("roc_auc"), pr_summary.rename("average_precision")],
    axis=1,
).style.format("{:.3f}")

	roc_auc	average_precision
label
short	0.969	0.826
hold	0.825	0.555
long	0.937	0.953
micro	0.932	0.877
macro	0.910	0.778

ROC curve

Each solid line treats one SPY regime as positive and the other two as negative. AUC measures how well the AAPL model ranks that regime across every decision threshold.

q.plot.roc(
    y_true,
    y_score,
    title="AAPL model of next-session SPY regime: ROC",
).show()

Precision-recall curve

Average precision summarizes precision across recall levels. This view makes performance on the less frequent short and hold regimes easier to judge than ROC alone.

q.plot.precision_recall(
    y_true,
    y_score,
    title="AAPL model of next-session SPY regime: precision-recall",
).show()