API Reference¶

def Problem(  # noqa: N802
    X: np.ndarray,  # noqa: N803
    target: np.ndarray | None = None,
    A: np.ndarray | None = None,  # noqa: N803
    b: np.ndarray | None = None,
    C: np.ndarray | None = None,  # noqa: N803
    d: np.ndarray | None = None,
    alpha: float = 0.0,
    rho: float = 0.0,
    mu: np.ndarray | None = None,
    target_lr: tuple[float, np.ndarray, np.ndarray] | None = None,
    pcg_lr: tuple[float, np.ndarray, np.ndarray] | None = None,
) -> _MinVarProblem | _Problem:
    """Create a portfolio optimisation problem.

    Returns a :class:`_MinVarProblem` (shrinking active-set) when no custom
    constraints are supplied, or a :class:`_Problem` (growing active-set) when
    any of ``A``, ``b``, ``C``, ``d`` are provided.

    Args:
        X:      Returns matrix of shape ``(T, N)``.
        target: Optional ``(N, N)`` regularisation matrix; when supplied the
                shrinkage term ``alpha * ||target @ w||^2`` is added to the
                objective.  ``None`` disables shrinkage entirely.
        A:      Equality constraint matrix ``(N, m)``: ``A^T w = b``.
        b:      Equality RHS ``(m,)``.
        C:      Inequality constraint matrix ``(N, p)``: ``C^T w <= d``.
        d:      Inequality RHS ``(p,)``.
        alpha:     Shrinkage intensity; only active when ``target`` is provided.
        rho:       Return tilt strength (Markowitz mean-variance).
        mu:        Expected returns vector ``(N,)``; required when ``rho != 0``.
        target_lr: Low-rank factored target ``(bar_lam, U_k, delta_k)`` for
                   RMT eigenvalue-cleaning; replaces ``target`` in the CG matvec.
        pcg_lr:    RMT preconditioner ``(bar_lam, U_k, delta_k)`` for
                   ``solve_pcg``; ignored unless PCG is invoked.

    Returns:
        A solver instance with ``solve_kkt()``, ``solve_minres()``,
        ``solve_cg()``, and ``solve_cvxpy()`` methods, each returning
        ``(w, n_iters)``.

    Examples:
        >>> import numpy as np
        >>> X = np.random.default_rng(42).standard_normal((500, 20))
        >>> w, _ = Problem(X).solve_kkt()
        >>> float(round(w.sum(), 8))
        1.0
        >>> bool((w >= 0).all())
        True
    """
    if A is None and b is None and C is None and d is None:
        return _MinVarProblem(X, target=target, alpha=alpha, rho=rho, mu=mu, target_lr=target_lr, pcg_lr=pcg_lr)

    # number of assets
    n = X.shape[1]

    A = A if A is not None else np.ones((n, 0))  # noqa: N806
    b = b if b is not None else np.ones(1)
    C = C if C is not None else -np.eye(n)  # noqa: N806
    d = d if d is not None else np.zeros(n)

    return _Problem(X, target=target, A=A, b=b, C=C, d=d, alpha=alpha, rho=rho, mu=mu)

def simulate_equity_returns(
    n: int,
    T: int,  # noqa: N803
    *,
    k: int | None = None,
    rng: np.random.Generator | int | None = None,
) -> np.ndarray:
    """Simulate a TxN demeaned equity return matrix with latent factor structure.

    Returns are generated from the model

        X = F @ B.T + E

    where F (TxK) are factor returns, B (NxK) are factor loadings, and E (TxN)
    is idiosyncratic noise.  The first factor is a market factor with universally
    positive loadings and high variance; the remaining k-1 factors are
    style/industry factors with sparse loadings.  This produces a covariance
    spectrum qualitatively similar to equity universes: a dominant market
    eigenvalue, a handful of secondary factor eigenvalues, and a long tail of
    near-equal idiosyncratic eigenvalues.

    Parameters
    ----------
    n:
        Number of assets.
    T:
        Number of time periods (trading days).
    k:
        Number of latent factors.  Defaults to ``max(3, n // 10)``.
    rng:
        Random state — a :class:`numpy.random.Generator`, an integer seed,
        or ``None`` (non-reproducible).

    Returns:
    -------
    X : ndarray of shape (T, n)
        Demeaned return matrix.  Each column has zero mean.

    Examples:
    --------
    >>> X = simulate_equity_returns(100, 200, k=5, rng=0)
    >>> X.shape
    (200, 100)
    >>> bool(abs(X.mean(axis=0)).max() < 1e-14)
    True
    """
    rng = np.random.default_rng(rng)
    if k is None:
        k = max(3, n // 10)

    # Factor volatilities (daily): market ~1 %, style factors ~0.5 %
    factor_vols = np.concatenate([[0.01], np.full(k - 1, 0.005)])

    # Factor returns: T x k
    F = rng.standard_normal((T, k)) * factor_vols  # noqa: N806

    # Factor loadings: n x k
    # Market: all assets have positive exposure in [0.4, 0.8]
    # Style:  sparse (~50 % non-zero), drawn from N(0, 0.2)
    B = np.zeros((n, k))  # noqa: N806
    B[:, 0] = rng.uniform(0.4, 0.8, size=n)
    for j in range(1, k):
        mask = rng.random(n) < 0.5
        B[mask, j] = rng.standard_normal(int(mask.sum())) * 0.2

    # Idiosyncratic volatility: uniform in [0.5 %, 1.5 %] per asset
    idio_vols = rng.uniform(0.005, 0.015, size=n)
    E = rng.standard_normal((T, n)) * idio_vols  # noqa: N806

    X: np.ndarray = F @ B.T + E  # noqa: N806
    X -= X.mean(axis=0)  # noqa: N806
    return X