arXiv CSopen access

Data Driven Block Replacement Scheduling

Data Driven Block Replacement Scheduling Aniruddhan Ganesaraman

Vidyadhar Kulkarni

aniruddhan [email protected]

[email protected]

Department of Statistics and Operations Research University of North Carolina at Chapel Hill Submitted to Operations Research June 27, 2026

Abstract.

arXiv:2607.15229v1 [cs.LG] 16 Jul 2026

We develop data-driven algorithms for maintaining N independent identical machines under a block replacement policy, in which each machine is replaced upon failure and all machines are jointly replaced at regular intervals of length k. The goal is to learn the cost-minimizing interval k ∗ from operational data when the lifetime distribution is unknown. At each decision epoch, the operator selects k ∈ {1, 2, . . . , K}, observes the resulting failure history (a mixture of complete and right-censored lifetimes) and incurs a per-unit-time cost governed by the renewal function. We formulate this as a stochastic multi-armed bandit and propose Hoeffding- and Bernsteinbased lower-confidence-bound algorithms achieving O(K log T ) regret, matching the Lai–Robbins lower bound. Exploiting a nested observation property unique to block replacement, correlated variants attain O((K − k ∗ ) log T ) regret and require only O(1) direct pulls of suboptimal arms k < k ∗ . A complementary Kaplan–Meier renewal algorithm estimates the lifetime distribution nonparametrically from censored data, achieving almost-sure policy consistency and empirically near-zero incremental regret at long horizons. We additionally analyze two average-cost MDPs: a time-elapsed formulation establishing that block replacement is optimal within its policy class for any lifetime distribution, and an age-vector formulation proving a monotone threshold structure under increasing failure rate distributions and providing a gold-standard cost benchmark. Numerical experiments confirm the theoretical ordering and reveal structural cost gaps between optimal block and age-dependent replacement. Code for replication is publicly available at this Github repository. Keywords. Block replacement; Multi-armed bandits; Censored lifetime data; Data-driven policy.

1.

Introduction

Many organizations face the problem of managing a large number of identical physical assets, all of which must remain functional at all times. For example, all the hard drives in a data center must remain operational to ensure continuous service availability. In aviation, every turbine blade in a commercial aircraft engine must remain in working condition for safe flight operations. Along 1

a stretch of railroad track, every fastening clip must remain secure to maintain safe train passage, and in ground vehicles, every brake pad must remain functional to meet safety standards. One simple way to achieve this is to replace each item when it fails. This involves dispatching a repair crew to the site and replacing the failed item(s). This is called an unplanned replacement. It can be expensive due to the fixed cost of mobilizing a crew for a single failure. One way to mitigate this is to replace some of the unfailed items when a crew visits to repair a failed one. However, this involves maintaining a list of which items to replace on each visit, and keeping such a list may be impossible or too costly and cumbersome. A simpler strategy is to replace all the items at the site at regular intervals, regardless of whether they have failed. This is called a block replacement policy. Thus each turbine blade in an aircraft engine is replaced when it fails, but every so often, in an engine overhaul (say every three years), all blades are replaced, whether failed or not. Likewise, each railroad fastener is replaced when found defective during inspection, but every so often (say every five years) all fasteners along an entire track section are replaced wholesale. The main advantage of a block replacement policy is that it is easy to implement, since it is described by a single parameter, the inter-replacement interval k. A block replacement policy with parameter k (denoted BR(k)) replaces each item upon failure and schedules a full replacement of all items every k time units. Clearly, the question is how to determine an optimal k, denoted by k ∗ . To determine this, we need a cost model that quantifies the costs involved. The simplest model is to assume that each failure (which leads to an unplanned replacement) costs cf dollars, while each block replacement costs cb dollars. Note that a smaller k reduces the cost of unplanned replacements, but increases the cost of block replacement. A larger k has the opposite effect. This implies that there is an optimal inter-replacement interval k ∗ that minimizes the long-run cost per unit time (cost rate). It is clear that the cost rate of the policy BR(k) depends on the lifetimes of the items. We consider the simplest model: the lifetimes of all items are independent and identically distributed (iid) random variables with a common cumulative distribution function (cdf) F . If F is known, it is easy to compute c(k), the long-run cost rate of BR(k), using renewal reward processes. See Section 3 for details. The problem of finding an optimal k then reduces to a one-dimensional optimization problem of minimizing c(k) over the allowed values of k. One can find an optimal k ∗ either analytically or numerically. This raises two fundamental questions. The first is how frequently to perform block replacements, i.e., how to determine k ∗ . The second, and subtler, question is what information should be used when making the replacement decision. At one extreme, the operator tracks only the time elapsed since the last block replacement, which is a single scalar. At the other extreme, the operator tracks the current age of every individual item, which is an N -dimensional vector. Time-elapsed policy. We formulate a Markov decision process (MDP) in which the state is the time n elapsed since the last block replacement, and the decision at each period is whether to perform a block replacement or to continue for one more period, assuming that F is known. We show that an optimal policy is a threshold policy at ks = k ∗ (Section 3.3). Crucially, this MDP has 2

a state space of cardinality K + 1 regardless of fleet size N , and directly grounds the data-driven learning algorithms described below. Age-vector policy. We also formulate an MDP in which the state is the current age of each item and the decision at each period is whether to replace all items now or to continue for one more period (in which case individual failures are replaced instantaneously), assuming F is known. We show that the optimal policy has a monotone threshold structure under increasing failure rate (IFR) distributions (Theorem 3.4), and that the MDP optimum can outperform c(k ∗ ) by as much as 50% (Section 5). Crucially, the state space grows as (K + 1)N , making this MDP computationally intractable for large fleets; it nonetheless serves as a gold-standard benchmark for ∗ between block and evaluating BR(k ∗ ) and quantifying the structural cost gap δstruct = c(k ∗ ) − gage

age-dependent replacement. We now address another important issue that arises during implementation. As mentioned before, the value of k ∗ depends on the lifetime distribution F . In practice, however, F is often unknown. We assume that we can keep records of when the unplanned and block replacements have occurred so far. We begin with N new items at time zero and schedule the first block replacement at time k1 , chosen in some reasonable fashion. At the n-th block replacement time, we use the data collected so far to determine the next inter-replacement time kn+1 . Such a policy is called a data-driven block replacement policy. We study several such policies in this paper. These data-driven policies can be broadly classified as one-step or two-step policies. One-step policies use the observed data to directly estimate c(k) and then choose the next k as a minimizer of this estimated c(k). Two-step policies use the current data to estimate the lifetime distribution F , then use this estimated F to compute c(k) and select the next k as its minimizer. We study both types. To keep the analysis tractable, we assume lifetimes are positive integer-valued random variables. (i) One-Step Policy: Bandit Algorithms. We assume there are a finite number of choices for k, say {1, 2, . . . , K}. We think of each choice as a bandit arm and follow the lower-confidencebound (LCB) algorithm to select a new arm at each decision point. Suppose we follow this policy for T decision points. We show that such a policy achieves O((K − 1) log T ) cumulative regret (Section 4.1), matching the Lai and Robbins [1985] lower bound. We improve this policy by exploiting a special feature of the problem: choosing interval k simultaneously reveals the costs (that would have been incurred) of all shorter intervals j ≤ k. This leads to bandit algorithms with correlated arms (Section 4.1.2) that require only O(1) direct pulls of suboptimal arms k < k ∗ , improving the regret bound from O((K − 1) log T ) to O((K − k ∗ ) log T ). (ii) Two-Step Policy: Kaplan–Meier Estimation.

At each block replacement we have

N machines that are replaced before they fail, producing N right-censored observations. We also observe a number of failures since the last block replacement, each yielding an uncensored observation. We use the Kaplan–Meier procedure Kaplan and Meier [1958] to estimate F from this 3

mix of censored and uncensored observations (Section 4.2). We use the estimated F to estimate c(k) and select b k ∗ as its minimizer. We prove that the estimated cost converges almost surely to c(k ∗ ) and the estimated optimal policy converges to k ∗ (Theorem 4.6). Numerically, this two-step procedure consistently outperforms the one-step bandit algorithms, albeit at higher computational cost per cycle. We view this work as a contribution to the reliability and maintenance literature that employs sequential learning tools, rather than a methodological contribution to bandit or learning theory in its own right. This paper is organized as follows. Section 2 reviews related work. Section 3 presents the model, regret framework, and both the MDP formulations. Section 4 develops the bandit and Kaplan–Meier algorithms. Section 5 reports numerical experiments. Section 6 states conclusions and directions for future research. Pseudocodes for all seven algorithms are collected in the Online Appendix.

2.

Literature Review

Multiple streams of literature are relevant to our work. We discuss them systematically below. We focus on the streams most directly relevant to our contributions. Classical block replacement.

The analytical foundations of block replacement were laid by

Barlow and Hunter [1960], who derived optimal replacement intervals under parametric lifetime distributions. Berg [1976] extended the framework to systems with increasing running costs. The close relationship between block replacement costs and the renewal function, including the subtlety that IFR does not imply convexity of the renewal function, is examined in detail by Shaked and Zhu [1992]. For multi-unit systems, Sheu and Jhang [1996] provide explicit cost formulae for group replacement with minimal repair. Cho and Parlar [1991] survey multi-unit maintenance models, and Dekker et al. [1997] focus specifically on economic dependence, the cost structure (simultaneous replacement is cheaper per machine than individual replacement) that motivates block replacement in this paper. MDP formulations and structural results.

Kurt and Kharoufeh [2010] established threshold

policies for a single machine with Markovian degradation, showing that the optimal replacement region is characterised by a state threshold. Andersen et al. [2022] compared MDP algorithms for multi-component replacement problems numerically, illustrating the computational challenges that motivate scalable data-driven alternatives. The theoretical foundations for average-cost MDPs used in Section 3.4, including Bellman optimality and convergence of value iteration, are developed in Puterman [1994] and Ross [1996]. Adaptive and data-driven maintenance.

The problem of learning maintenance policies from

operational data has grown substantially in importance. 4

Corman et al. [2017] used Kaplan–

Meier estimates of failure time data to calibrate preventive maintenance schedules for railway braking systems, and Deprez et al. [2021] combined survival regression with renewal theory to price maintenance contracts. Most directly related to the present paper is Ge and Kulkarni [2024], who formulate single-machine maintenance as a renewal-reward bandit and develop data-driven replacement scheduling algorithms. Our work extends this framework to block replacement policies, introduces variance-adaptive and correlated-arm LCB variants, demonstrates that our algorithms match the Lai–Robbins lower bound, and adds the Kaplan–Meier renewal approach. Multi-armed bandit algorithms.

Lai and Robbins [1985] established the fundamental Ω(log T )

lower bound on regret for any consistent algorithm; Auer et al. [2002] proved that UCB achieves O(log T ) regret for bounded rewards. Audibert et al. [2009] introduced variance-adaptive algorithms based on the empirical Bernstein inequality, yielding tighter instance-dependent constants when arm reward variances differ; we adapt these for our Bernstein-based LCB algorithms. Side observations and correlated arm structures were studied by Mannor and Shamir [2011] and Caron et al. [2012]. The cascade bandit of Kveton et al. [2015], where a single play reveals a cascade of outcomes, is the closest structural analogue to our nested observation property. A comprehensive survey of bandit algorithms is given in Bubeck and Cesa-Bianchi [2012]. Nonparametric renewal estimation under censoring.

The Kaplan–Meier estimator Kaplan

and Meier [1958] is the nonparametric maximum likelihood estimator for survival functions under right censoring. Vardi [1982] first studied nonparametric estimation in renewal processes from multiple independent streams; Grübel and Pitts [1993] constructed an empirical renewal function √ estimator for complete data and proved n-consistency. Under censoring, Baxter and Li [1995] showed that the KM plug-in estimator of the renewal function is consistent under mild conditions. Our KM algorithm in Section 4.2 builds directly on these results. Our paper adapts and extends these existing methodologies to the specific problem of deriving optimal block replacement policies. To the best of our knowledge, this work is novel and has not appeared in literature before.

3.

Model Specification

3.1

Model Description

Consider a system of N identical machines operating independently of each other. Each machine j ∈ {1,  ·· · , N } has independent, identically distributed (possibly unbounded) integer-valued lifetimes Lji , where Lji represents the ith lifetime of machine j. Let L denote a generic random variable i≥1

with the same common distribution with probability mass function f : N+ → [0, 1] and cumulative t P distribution function F : N+ → [0, 1], where F (t) = f (s) = P(L ≤ t) for t = 1, 2, · · · . s=1

We impose the following structural assumptions about the lifetimes: 5

(i) The lifetime sequences {Lji : i ≥ 1} for different machines j = 1, . . . , N are mutually independent and identically distributed taking values on {1, 2, 3, · · · }. This is without essential loss of generality, as continuous lifetimes can be discretized to match any reporting resolution. (ii) The lifetime distribution has finite mean µ := E[L] =

∞ P

t · f (t) < ∞.

t=1

The system operates under a block replacement policy that combines scheduled block replacements with individual replacements immediately upon failure. In block replacements, all N machines are replaced simultaneously at a fixed cost cb > 0, regardless of their operational status. When a machine fails between block replacements, it is immediately replaced at cost cf > 0. If a machine fails exactly at a scheduled block replacement time, such a failure is recorded as an individual replacement at cost cf and the block replacement at cost cb still proceeds as scheduled. All replacements (block and individual) are assumed to occur instantaneously. Upon replacement, a machine begins a new lifetime with distribution F . Machines can be replaced at any discrete time instant regardless of its operational status. We assume cf > cb /N , which ensures that individual replacements are more expensive per machine than block replacements, thereby incentivizing the use of block replacement policies. This cost structure incentivizes proactive block replacement but must be balanced against unnecessary replacement of functional machines. Let τ0 = 0 and let τt denote the tth block replacement time for t ≥ 1. Define the inter-block replacement interval as kt := τt − τt−1 for t ≥ 1. The tth cycle of the system spans the time instants {τt−1 + 1, . . . , τt } for t ≥ 1. For each machine j ∈ {1, . . . , N }, let Ytj denote the (random) number of failures (and hence individual replacements) of machine j during the tth cycle. The total cost N P incurred during the t-th cycle is ψt = cb + cf · Ytj . For t ≥ 1, the cost per unit time incurred during the tth cycle is given by γt = ψktt =

j=1 N P cb +cf Ytj j=1

kt

.

Let (N (t))t≥0 denote the renewal process with interarrival distribution F , where N (t) represents the number of renewals (failures) in the interval [0, t]. The associated renewal function is defined as M (t) := E[N (t)]. As shown in Kulkarni [2009], M (t) satisfies the renewal equation: M (t) = F (t) +

t X

f (s)M (t − s),

(1)

s=1

with initial condition M (0) = 0. A block replacement policy with parameter k ∈ N+ , denoted BR(k), is a stationary policy under which all machines are replaced every k time units, regardless of their operational status. Between block replacements, machines are replaced individually upon failure. Under BR(k), the system operates in cycles of deterministic length k. By the renewal reward theorem (see Ross [1996]), the in one cycle] = long-run average cost per unit time is given by c(k) = E[Cost E[Cycle length]

6

cb +cf ·N ·M (k) , because the k

expected number of failures for a single machine over the interval [0, k] is M (k), and we have N independent machines. The set of feasible replacement intervals is [K] = {1, 2, . . . , K}, where K ∈ N+ is a known upper bound on block replacement intervals. This bound may arise from operational constraints like maximum time between scheduled maintenance. Further we assume that f (k) > 0 for all k ∈ [K], that is, the machines can fail after reaching age k for any k ∈ [K]. The optimal block replacement interval is defined as k ∗ := arg mink∈[K] c(k). The corresponding optimal long-run average cost is c∗ := c(k ∗ ). For simplicity, we assume that k ∗ is unique, that is, c(k ∗ ) < c(k) for all k ̸= k ∗ . This assumption does not affect the analytic interpretability of the results that follow. Remark 3.1. The cost function c(·) need not be convex in general (see Section 5). Convexity of the renewal function M (·) is sufficient to ensure convexity of c(·), but unfortunately, the convexity of the renewal function M (·) is not characterized by simple properties of the lifetime distribution F . In particular, increasing failure rate (IFR) of L does not imply convexity of M (·). See Shaked and Zhu [1992] for a detailed discussion on this topic.

3.2

Regret Framework for Online Learning

In practical applications, the lifetime distribution F (and hence the renewal function M and optimal policy k ∗ ) is typically unknown. We must learn the optimal policy from observed data while the system operates. The unknown F motivates an online learning framework. Consider an algorithm that produces a sequence of replacement interval decisions {kt }t≥1 , where kt ∈ [K] is the block replacement interval chosen for the tth cycle. For k ∈ [K], let ∆k := c(k)−c(k ∗ ) denote the suboptimality gap. The performance of the algorithm is evaluated by the cumulative regret over a horizon of T cycles, defined as ∗

Reg(T ) := E[Z(T )] − T · c(k ) =

T X

E[∆kt ],

(2)

t=1

where Z(T ) :=

T P

γt is the sum of per-unit-time costs incurred during the first T cycles. This

t=1

measures the expected excess per-unit-time cost compared to an oracle that always uses the optimal policy k ∗ . We also introduce an alternative formulation that is more interpretable in terms of total cost incurred. Define T h i X ∗ g e Reg(T ) := E Z(T ) − c(k )S(T ) = kt E[∆kt ],

(3)

t=1

e )= where Z(T

T P

ψt is the total cost incurred, and S(T ) =

t=1

T P t=1

7

kt is the total elapsed time, over the

first T cycles. Since 1 ≤ kt ≤ K, the two versions of regret are related by g ) ≤ K · Reg(T ), Reg(T ) ≤ Reg(T

(4)

so both grow at the same asymptotic rate in T . We primarily use Reg(T ) for our theoretical analysis, while both versions are compared in our numerical studies. Regret is analysed by its asymptotic growth rate. We write Reg(T ) = O(f (T )) if there exist positive constants C and T0 such that for all T ≥ T0 , we have Reg(T ) ≤ Cf (T ). Similarly, for two functions f, g : N → R+ , we say that g is Ω(f ) if there exist constants c > 0 and n0 ∈ N such that g(n) ≥ c · f (n) for all n ≥ n0 . Before presenting our algorithms, we recall fundamental limits on achievable regret due to Lai and Robbins [1985]. Let Pk be the distribution of the observed cost under policy BR(k). With DKL (· ||·) denoting the Kullback-Leibler divergence between two distributions, Lai & Robbins Lai and Robbins [1985] show that for any “consistent” algorithm, X E[Reg(T )] ∆k ≥ . T →∞ log T DKL (Pk ∥Pk∗ )

lim inf

(5)

k∈[K]: k̸=k∗

g ) = Ω(log T ). The key insight is that This says Reg(T ) = Ω(log T ) and hence by (4), Reg(T distinguishing policy k from k ∗ requires Ω(log T /KL(Pk ∥Pk∗ )) observations, so the best achievable regret grows as O(log T ). Z(T )] ∗ e Remark 3.2. For any sequence {kt }t≥1 of decisions, we have (i) E[ E[S(T )] ≥ c(k ) because E[Z(T )] = T T P P )] E[kt c(kt )] ≥ c(k ∗ ) E[S(T )] and (ii) E[Z(T ≥ c(k ∗ ) because E[Z(T )] = c(kt ) ≥ T · c(k ∗ ). T e

t=1

3.3

t=1

Time-elapsed MDP Formulation

This Markov Decision Process (MDP) uses the time elapsed since the last block replacement as the state variable. This yields a state space of cardinality K + 1, independent of fleet size N . We assume the lifetime distribution F is known. Let k ∈ [K] denote the number of time periods elapsed in the current block replacement cycle, and let k = 0 denote the state immediately following a block replacement. Define the renewal mass function m(k) := M (k) − M (k − 1) for k ≥ 1, where M (k) = E[N (k)] is the renewal function of Section 3. The quantity m(k) equals the expected number of failures per machine during the k th period of a renewal process in which each machine is replaced immediately upon failure; the expected total number of failures across the fleet of N independent machines during period k is N m(k). The action space is A = {0, 1}; at each state k ∈ [K], the operator selects one of two actions: 8

continue (a = 0), incurring expected cost cf N m(k) from individual failure replacements and transitioning to state k + 1; or block replace (a = 1), incurring cost cb and transitioning to the reset state k = 0, from which the next cycle begins immediately at state k = 1. Block replacement is mandatory at k = K since K is the upper bound on feasible replacement intervals. The expected one-step costs are c(0, k) = cf N m(k) and c(1, k) = cb . We adopt the undiscounted average cost criterion for this infinite-horizon problem.

Under a

stationary deterministic policy π : [K] → {0, 1}, the long-run average cost per unit time starting from state k is

"T −1 # 1 π X l (k) = lim sup E c(π(kt ), kt ) k0 = k . T →∞ T π

(6)

t=0

The optimal average cost is l∗ = inf

π

lπ (k). As shown in Ross [1996], the optimal long-run average

cost l∗ and the relative value function v : {0, . . . , K} → R, normalized so that v(0) = 0, satisfy the average-cost Bellman optimality equation for k ∈ [K] given by
l∗ + v(k) = min cf N m(k) + v(k + 1), cb + v(0) .

(7)

An optimal policy π ∗ achieves the minimum in (7) at every state. Algorithm 1 estimates l∗ and v by iterating the Bellman operator. We set k = 0 as the reference state, so w(n) (0) = 0 for all n, and monitor convergence via the span seminorm sp(d(n) ) := maxk d(n) (k) − mink d(n) (k) where d(n) (k) := w(n+1) (k) − w(n) (k). Theorem 3.1 (Convergence of Value Iteration). The time-elapsed MDP satisfies the unichain condition: state 0 is accessible from every state k ∈ [K] under the block replacement action. Consequently, the value iteration sequence {w(n) } generated by Algorithm 1 converges, and upon termination with sp(d(n) ) < ε, the following hold: (i) the midpoint estimate satisfies |b l∗ − l∗ | ≤ ε/2; (ii) the computed policy π b∗ is ε/2-optimal. Proof. Proof. This is a standard result for unichain average-cost MDPs; see Puterman [1994], Theorem 8.5.6.

The unichain condition holds because the block replacement action (a = 1)

transitions any state directly to state 0.

■

A stationary deterministic policy π is characterized by its replacement set Sπ := {k : π(k) = 1}. Since the state resets to 0 upon any block replacement and then advances deterministically through 1, 2, . . ., each cycle under π runs from state 0 until the system first enters Sπ . States beyond min(Sπ ) are never visited within a single cycle: once the system reaches min(Sπ ), it replaces and restarts. Consequently, if min(Sπ ) = k, the cycles are i.i.d. with deterministic length k and expected cost cb + cf N M (k), so the long-run average cost is lπ = c(k) by the renewal reward theorem. In particular, lπ depends on π only through min(Sπ ). Theorem 3.2. For every lifetime distribution F with E[L] < ∞, the following hold for the timeelapsed MDP: 9

(i) The optimal average cost is l∗ = c(k ∗ ). (ii) The threshold policy πk∗ (k) = 1[k ≥ k ∗ ] is an optimal stationary policy.  (iii) The threshold ks := min k ∈ [K] : cf N m(k) + v(k + 1) > cb + v(0) = k ∗ . Proof. Proof. Note that under any stationary deterministic policy π, the system resets to state 0 upon every block replacement and then advances through states 1, 2, · · · deterministically until it first enters Sπ . Each cycle therefore has deterministic length min(Sπ ) and expected cost cb + cf · N M (min(Sπ )), since the expected number of individual failures per machine over an interval of length min(Sπ ) is M (min(Sπ )) by definition of the renewal function. The renewal reward theorem (Ross 1996) then gives lπ =

cb +cf ·N M (min(Sπ )) min(Sπ )

= c(min(Sπ )). Since lπ = c(min(Sπ )) for every

stationary deterministic policy π, optimizing over all such policies reduces to minimizing c(k) over k ∈ [K]. The infimum is therefore c(k ∗ ), achieved by any policy with min(Sπ ) = k ∗ ; the threshold policy πk∗ is one such policy. By Theorem 8.5.4 of Puterman [1994], there exists an optimal stationary deterministic policy, so l∗ = c(k ∗ ) and πk∗ is optimal. For part (iii): under πk∗ with l∗ = c(k ∗ ), the Bellman equation (7) gives cf N m(k) + v(k + 1) ≤ cb + v(0) for all k < k ∗ (continue is optimal) and cb + v(0) ≤ cf N m(k ∗ ) + v(k ∗ + 1) (block replace is optimal at k ∗ ), so ks = k ∗ . ■ While the time-elapsed MDP achieves the optimal cost c(k ∗ ) within the class of block replacement policies, it discards information about the current ages of individual machines. A richer formulation that tracks the full age vector of the fleet can exploit this information to achieve a strictly lower average cost under IFR lifetime distributions. We develop this formulation next.

3.4

Age-vector MDP Formulation

This MDP approach models the system state explicitly by tracking the ages of all machines, enabling the computation of optimal policies. This framework allows for state-dependent policies that can exploit information about the current operational ages of individual machines. Throughout this subsection, we assume that the lifetime distribution F is known. This MDP formulation serves as a benchmark to evaluate the performance gap between block replacement policies and fully state-aware policies. In Section 5, we numerically compare the costs achieved by the proposed bandit (Section 4.1) and Kaplan - Meier (Section 4.2) algorithms (which learn without knowing F ) with this benchmark. j Let the system state at the beginning of time step n be sn = (s1n , . . . , sN n ) ∈ S, where sn ∈

{0, 1, . . . , K} represents the current age (time since last replacement) of machine j. The state space is S = {0, 1, . . . , K}N , with cardinality |S| = (K + 1)N . Let S(t) = 1 − F (t) denote f (t) the survival function, h(t) := S(t−1) = P (L = t | L ≥ t) the discrete hazard function, and

h̄(t) := 1 − h(t) = S(t)/S(t − 1) be the survival ratio for t ≥ 1. We impose the boundary condition h(K + 1) = 1 so that h̄(K + 1) = 0. This means that a machine which has survived K units of 10

time since its last replacement fails with certainty at the next step; age K + 1 is therefore never reached in the state space. The action space is A = {0, 1}, where action a = 0 denotes continuation (individual machines replaced on failure) and a = 1 denotes block replacement of all machines, transitioning the system to 0 := (0, . . . , 0). Under action a = 0, each machine j at age sj fails over the next period independently with probability h(sj + 1) and survives with probability h̄(sj + 1). A failed machine is immediately replaced (age resets to 0); a surviving machine ages by one unit to sj + 1. Since h(K + 1) = 1, if machine j survives, then sj ≤ K − 1. The transition therefore never produces an age exceeding K. Since machines fail independently, the set of machines that fail over one time step is a random subset A ⊆ [N ]. For each subset A ⊆ [N ] := {1, 2, · · · , N } of machines that fail, the resulting next state s[A] ∈ S is defined by  0 s[A],j := sj + 1

j ∈ A,

(8)

j∈ / A.

Q For A ⊆ [N ], the one-step transition probabilities under action a = 0 are P0 s, s[A] = j∈A h(sj +1)· Q j ′ ′ [A] for any A ⊆ [N ]. Under action a = 1 j ∈A / h̄(s +1), with P0 (s, s ) = 0 if s is not of the form s (block replacement), the transition is deterministic: P1 (s, s′ ) = 1s′ =0 . The expected one-step costs are given by c(0, s) = cf

N P

h(sj + 1) and c(1, s) = cb . The cost under

j=1

a = 0 follows because machine j fails independently with probability h(sj + 1), and by linearity of N P expectations, the expected number of failures across N machines is h(sj + 1). j=1

We adopt the undiscounted average cost criterion for this infinite horizon problem. Under a stationary deterministic policy π : S → A, the long-run average cost per unit time starting from state s is defined as

"T −1 # X 1 g π (s) = lim sup Eπ c(π(st ), st ) s0 = s . T →∞ T

(9)

t=0

where the expectation is taken over state transitions induced by policy π. ∗

The optimal average cost is g ∗ = inf π g π (s), and a policy π ∗ is optimal if g π (s) = g ∗ for all s ∈ S. As proved in Ross [1996], the optimal average cost g ∗ and relative value function w∗ : S → R satisfy the Bellman optimality equation for each s ∈ S, given by n o X g ∗ + w∗ (s) = min c(a, s) + Pa (s, s′ )w∗ (s′ ) . a∈A

(10)

s′ ∈S

The relative value function captures the advantage (or disadvantage) of starting in state s versus the reference state, in terms of accumulated cost before the process regenerates. An optimal policy π ∗ achieves the minimum in (10) at every state. Algorithm 2 describes value iteration to estimate g ∗ 11

and w∗ by iterating the Bellman operator. We set sref = 0 is the reference state (so w(k) (0) = 0 for all k). The convergence criterion is using the span seminorm sp(d(k) ) := maxs d(k) (s) − mins d(k) (s), where d(k) (s) := w(k+1) (s) − w(k) (s), which is standard for average cost MDPs. The reference state sref = 0 is arbitrary; any state can be used since the relative value function is unique only up to an additive constant.

Theorem 3.3 (Convergence of Value Iteration). The MDP described above satisfies the unichain condition, that is, there exists a state that is accessible from every other state under every stationary policy. Consequently, the value iteration sequence {w(k) } generated by Algorithm 2 converges, and upon termination with sp(d(k) ) < ε, the following hold: (i) the midpoint estimate satisfies |b g ∗ − g ∗ | ≤ ε/2; ∗

(ii) the computed policy π b∗ is ε/2-optimal, i.e., g πb − g ∗ ≤ ε/2. Proof. Proof. This is a standard result in MDP literature (for instance, see Puterman [1994]). The unichain condition holds because state 0 is accessible from any state s ∈ S because the f (k) > 0 for all k ∈ [K] and so the probability of all machines failing at a given time is positive.

■

Remark 3.3. The average-cost optimality equation (10) holds over all policies, including historydependent ones. The unichain property (Theorem 3.3) then guarantees that the infimum long run average cost is achieved by a stationary Markovian policy π ∗ (see Ross [1996]). The optimal policy for the MDP formulation possesses several structural properties which are summarized in the following theorem. Equip S with the componentwise partial order as s ⪰ s′ if sj ≥ s′j for all j ∈ [N ].

Theorem 3.4. Let π ∗ denote an optimal policy for the MDP. Suppose the lifetime distribution F has increasing failure rate (IFR), that is, h(t) is non-decreasing in t. Then the following are true: (a) Value Function Monotonicity: The optimal relative value function w∗ is non-decreasing in each coordinate, that is, w∗ (s′ ) ≥ w∗ (s) whenever s′ ⪰ s. (b) Policy Monotonicity: If π ∗ (s) = 1 and s′ ⪰ s, then π ∗ (s′ ) = 1. (c) Symmetry: The optimal relative value function and optimal policy are invariant under
permutation of machine indices, that is, for any permutation σ of [N ], w∗ s1 , . . . , sN =


w∗ sσ(1) , . . . , sσ(N ) and π ∗ s1 , . . . , sN = π ∗ sσ(1) , . . . , sσ(N ) . (d) Threshold Structure: For each fixed i ∈ [N ] and s−i ∈ {0, . . . , K}N −1 , there exists a threshold τi∗ (s−i ) ∈ {0, . . . , K} ∪ {∞} such that π ∗ (si , s−i ) = 1si ≥τi∗ (s−i ) , where we set τi∗ (s−i ) = ∞ if π ∗ (·, s−i ) ≡ 0. Moreover, τi∗ (s−i ) is non-increasing in each component of s−i . 12

The proof is in appendix A. Remark 3.4. Theorem 3.4 establishes that the optimal policy is completely characterised by the N threshold functions {τi∗ (s−i )}i∈[N ] , requiring N (K + 1)N −1 values in total, rather than a binary decision at each of the (K + 1)N states. By the symmetry of part (c), τi∗ (s−i ) depends only on the multiset of ages {s1 , . . . , sN }, so one need only store threshold values for states with
s1 ≤ s2 ≤ · · · ≤ sN , of which there are K+N , reducing storage by a factor of up to N !. N Remark 3.5. The time-elapsed MDP (Section 3.3) has a state space of cardinality K + 1, requires no assumption on F beyond E[L] < ∞, is tractable for all fleet sizes N , and establishes that block replacement is optimal within the class of time-elapsed policies, achieving cost g ∗ = c(k ∗ ). The age-vector MDP of Section 3.4 has a state space of cardinality (K + 1)N , requires IFR for its ∗ ≤ c(k ∗ ) structural results, is tractable only for small N , and achieves the lower optimal cost gage ∗ ≥ 0 by exploiting the joint distribution of machine ages. The structural gap δstruct = c(k ∗ ) − gage

quantifies the value of age information and is strictly positive under IFR for N ≥ 2.

4.

Data-Driven Algorithms

4.1

Multi-arm bandit based

We formulate the block replacement problem as a stochastic multi-armed bandit (MAB) with K arms, where each arm k ∈ [K] represents a block replacement interval of length k. At each decision epoch t ≥ 1, the decision maker selects a block replacement interval kt , observes the total number of failures Yt that occur across all N machines during that interval, incurring a total cost of cb + cf Yt   in that cycle. Let µk := E C (k) = cb + cf · N · M (k) = k · c(k) for k ∈ [K]. The goal is to identify k ∗ = arg mink∈[K] µkk , the optimal arm that incurs the minimum cost per unit time. The per-unit-time cost C (k) /k has the same range (width of support) b := cf N for every arm k ∈ [K], since the minimum is cb /k (no failures) and the maximum is cb /k + cf N (every machine fails every period). Throughout, log denotes the natural logarithm. 4.1.1

Independent Arms

In the independent arms setting, each arm k is treated as an independent learning problem. We develop lower confidence bound (LCB) algorithms based on Hoeffding’s and Bernstein’s concentration inequalities. 4.1.1.1

Hoeffding-based.

If k has been pulled nk (t) times up to decision epoch t, with observed nP (k) k (t) Ci (k) (k) costs C1 , . . . , Cnk (t) , the empirical mean cost associated to arm k ∈ [K] is b ck (t) = nk1(t) k . i=1

The Hoeffding-based lower confidence bound for arm k is s LCBIH ck (t) − b k (t) = b 13

2 log(t) . nk (t)

The algorithm selects the arm with the minimum LCB, that is kt = arg mink∈[K] LCBIH k (t − 1). Theorem 4.1. Under the Hoeffding-based algorithm (Algorithm 3) with independent arms, the expected cumulative regret satisfies  Reg(T ) = O 

X b2 log T

k̸=k∗

∆k

  = O ((K − 1) log T ) .

Proof. Proof. This is a standard result in MAB literature that follows from the arm-wise decomposition of regret given in Equation 2 and the fact that E[nk (T )] = O(log T ) for k ̸= k ∗ (see Lattimore and Szepesvári [2020]). 4.1.1.2

■

Bernstein-based.

Bernstein’s inequality provides tighter confidence bounds (smaller

exploration bonus) by exploiting variance information. For arm k ∈ [K] with nk (t) ≥ 2 samples, we compute the empirical variance as nk (t)

σ bk2 (t) =

!2

(k)

Ci k

X 1 nk (t) − 1 i=1

−b ck (t)

.

ck (t) − The Bernstein-based LCB is constructed as follows: when nk (t) = 1, define LCBIB k (t) = b p b 2 log t /nk (t); when nk (t) ≥ 2, set ck (t) − LCBIB k (t) = b

q 2.4 σ bk2 (t) log t /nk (t)

− 3.6b log t /nk (t)

(11)

The specific constants 2.4 and 3.6 are chosen to match the analysis carried out in Audibert et al. IH [2009]. The algorithm proceeds very similar to Algorithm 3, but uses LCBIB k (t) instead of LCBk (t).

Theorem 4.2. Under the Bernstein-based algorithm (Algorithm 4) with independent arms, the expected cumulative regret satisfies  Reg(T ) = O 

X b2 log T

k̸=k∗

∆k

  = O ((K − 1) log T ) .

Proof. Proof. See Theorem 4 of Audibert et al. [2009]. 4.1.2

■

Correlated Arms

A key structural property of the block replacement problem is the nested observation structure: when selecting a block replacement interval k and observing the failure history up to time k, we simultaneously obtain complete information about what would have happened under all smaller 14

intervals j < k. This “correlation” among arms can be exploited to improve learning efficiency. The following proposition formalizes this observation.

Proposition 4.1 (Nested Observations). When arm kis pulled, for each l ≤ k define C (l) :=  N n P j P cb + cf N (l) where N j (l) := max n ≥ 0 : Lji ≤ l . j=1

i=1

(a) C (l) is determined by the observed failure history of the pull of arm k; d

e (l) , where C e (l) is the cost from an independent direct pull of arm l. (b) C (l) = C Proof. Proof. (a) Since N j (l) ≤ N j (k) a.s., N j (l) is a deterministic function of (Lj1 , . . . , LjN j (k) ), which are observed when arm k is pulled. (b) The joint law of (Lj1 , . . . , LjN j (l) )j∈[N ] is identical whether these lifetimes are extracted from a pull of arm k ≥ l or from an independent direct pull e (l) are the same function of identically distributed inputs. of arm l, so C (l) and C ■ Let mk (t) denote the number of samples collected for arm k up to time t under the correlated arms scheme. When arm k is pulled at time t, we get mj (t) = mj (t − 1) + 1 samples for all j ≤ k upto time t. The number of samples for arm k equals the total number of times any arm j ≥ k has been K P nj (t). For any arm k, we have mk (t) ≥ nk (t), with equality only if no pulled, that is mk (t) = j=k

arm larger than k has ever been pulled. The sample counts satisfy m1 (t) ≥ m2 (t) ≥ · · · ≥ mK (t) for all t, because mk (t) − mk+1 (t) = nk (t) ≥ 0. Smaller arms always accumulate at least as many samples as larger arms.

Remark 4.1. Proposition 4.1 justifies the simultaneous update of empirical means b cj for all j ≤ kt upon pulling arm kt because it guarantees that these derived observations are equivalent to direct K P pulls of arm j, so the effective sample count for arm j is mj (t) = nl (t) ≥ nj (t). l=j

Remark 4.2. The regret contribution from arm k is ∆k · E[nk (T )], and not ∆k · E[mk (T )], because only direct pulls of arm k incur the suboptimality cost ∆k . The “free samples” contribute to learning the cost but not to regret, which is why the correlated algorithm achieves better regret bounds.

4.1.2.1

Hoeffding-based.

Using the nested observations structure, we update the sample

means for all arms k ≤ kt when arm kt is selected using the equation s (m) LCBCH ck (t) − b k (t) = b

15

2 log(t) , mk (t)

(m)

(m)

where b ck (t) is the empirical mean based on the mk (t) samples defined by b ck (t) = mk1(t)

mP k (t) i=1

(k)

Ci k

.

Theorem 4.3. Under the Hoeffding-based algorithm (Algorithm 5) with correlated arms, the expected cumulative regret satisfies  Reg(T ) = O 

X b2 log T k̸=k∗

∆k

  = O ((K − k ∗ ) log T ) .

The proof is in Appendix B. 4.1.2.2

Similar to the independent arms case, the Bernstein-based LCB r 2 log t (m) with correlated arms uses LCBCB ck (t) − b when mk (t) = 1; and when mk (t) ≥ 2, k (t) = b mk (t) we set v   u (m) 2 u 2.4 σ b (t) log t t k (m) (t) = b c (t) − LCBCB k k mk (t) 3.6b log t (12) − mk (t) where



Bernstein-based.

 mP k (t) (m) 2 1 σ bk (t) = mk (t)−1 i=1



(k)

Ci k

2 −b ck (t) is the sample variance based on mk (t) ≥ 2.

Theorem 4.4. Under the Bernstein-based algorithm (Algorithm 6) with correlated arms, the expected cumulative regret satisfies  Reg(T ) = O 

X b2 log T k̸=k∗

∆k

  = O ((K − k ∗ ) log T ) .

The proof is in Appendix C. Remark 4.3. Algorithms 5 and 6 can be initialized by pulling arm K, so that every arm satisfies mj ≥ 1 and the initialization phase is complete in one cycle. Remark 4.4 (Online variance computation). Algorithms 4 and 6 compute σ bk2 from three running nk nk P P scalars (nk , Sk , Qk ) via σ bk2 = (Qk − Sk2 /nk )/(nk − 1), where Sk = Ci and Qk = Ci2 . This i=1

i=1

requires O(1) storage and O(1) arithmetic per update.

Remark 4.5. The lower bound (Section 3.2) described in Lai and Robbins [1985] requires that any consistent algorithm pulls each suboptimal arm at least Ω(log T ) times. The independent 16

algorithms (Algorithms 3–4) achieve O(log T ) pulls for every suboptimal arm k ̸= k ∗ , matching this lower bound up to constants. The correlated algorithms improve upon this for arms k < k ∗ : these arms receive O(1) direct pulls because each pull of any larger arm j > k provides a free sample for arm k (Proposition 4.1).

4.2

Survival Analysis based

Under a block replacement scheduling for machine maintenance, when a block replacement occurs at time k, every machine that has not yet failed only provides a censored observation of its lifetime because we do not observe its actual failure time. With traditional approaches that simply collect failure times and build an empirical distribution from complete observations alone, machines with longer lifetimes are disproportionately censored, leading to a biased sample that skews toward shorter lifetimes. This bias can result in suboptimal policies that replace machines too frequently. To address this, we develop an approach that properly accounts for censored observations using survival analysis techniques, specifically the Kaplan-Meier estimator, combined with renewal theoretic cost computation. Consider the tth cycle with inter-block replacement interval kt = k. Machine j ∈ [N ] generates Ytj uncensored observations (ℓ, 1) corresponding to the complete lifetimes Lj1 , . . . , Lj j , and one censored Yt

observation (ℓ, 0) for the residual time from its last failure to the block replacement epoch, since every machine is replaced at time k regardless of its operational status. Aggregating over all N P j machines, cycle t contributes N j=1 (Yt + 1) observations. The accumulated dataset after T cycles T is DT = {(ℓi , δi )}ni=1 , where

nT =

T X N X

(Ynj + 1) = N T +

n=1 j=1

T X N X

Ynj ,

n=1 j=1

with expected size N (M (k) + 1) per cycle of length k.

Kaplan-Meier (KM) Estimation Given the dataset DT with censored observations, we estimate the lifetime distribution using the Kaplan-Meier (KM) estimator, which is the nonparametric maximum likelihood estimator for the survival function in the presence of right censoring. Let t1 < t2 < · · · < tnT denote the distinct observed times (both failures and censorings) in DT . At each time tj , let dj denote the number of failures at time tj and let nj be the number of machines that have not yet failed just before time tj .   Q d b The Kaplan-Meier estimator of the survival function is S(t) = j:tj ≤t 1 − njj for t ≥ 0, with b b and fb(t) = S(t b − 1) − S(t) b and compute M c via the renewal equation S(0) = 1. Set Fb(t) = 1 − S(t) 17

(1) setting c(t) = Fb(t) + M

t X

c(t − s) fb(s)M

s=1

c(0) = 0. Then estimate b for 1 ≤ t ≤ K and M c(k) = ∗ b k = arg min b c(k).

c(k) cb +cf ·N ·M . The estimated optimal policy is k

k∈[K]

Algorithm 7 has two hyperparameters: an exploration rate ε and a refit frequency R. A coldstart phase (t ≤ 5) ensures initial KM coverage; thereafter the algorithm uses a decaying schedule √ k ∗ otherwise. The parameter εt = ε/ t, exploring uniformly with probability εt and exploiting b R controls how often the KM estimator and renewal function are recomputed: smaller R yields faster adaptation to accumulating data at higher per-cycle computational cost; larger R amortizes computation but delays policy updates. We now establish the consistency and convergence properties of the survival-based renewal approach.

Theorem 4.5 (Consistency of Kaplan-Meier Estimator ). Let FbT denote the Kaplan-Meier estimator based on dataset DT with nT observations. Suppose the censoring mechanism satisfies the independent censoring assumption. Then a.s. sup |FbT (t) − F (t)| −−→ 0

as nT → ∞.

t∈[K]

The proof is in Appendix D. Theorem 4.6 (Policy and cost convergence). Let b kT∗ denote the estimated optimal policy based on dataset DT with nT observations, and let k ∗ denote the true optimal policy. Under the conditions of Theorem 4.5, as nT → ∞, P(b k ∗ = k ∗ ) → 1. Furthermore, the estimated cost converges, that is a.s. |b c(b kT∗ ) − c(k ∗ )| −−→ 0.

T

The proof is in Appendix E.

4.3

Computational Complexity Analysis

The four bandit algorithms (Algorithms 3 – 6) share the same asymptotic computational complexity. Each maintains two or three running scalars per arm: (nk , Sk ) for Hoeffding and (nk , Sk , Qk ) for Bernstein, requiring O(K) space. Per decision epoch, computing all confidence bounds and selecting the minimum costs O(K); the correlated update (accumulating statistics for all j ≤ kt ) adds O(K) to time complexity per epoch. Thus total time complexity is O(KT ) for all four algorithms. The Kaplan-Meier algorithm maintains vectors c(·), F (·), and M (·) of length K, requiring O(K) c(t) via the renewal equation requires space. Updating c(·) and F (·) costs O(K), but computing M an O(K 2 ) convolution. The per-epoch cost is therefore O(K 2 ), giving O(K 2 T ) total. Despite 18

the higher cost, the KM algorithm correctly incorporates all censored observations and in practice identifies k ∗ with far fewer cycles than the bandit methods. Value iteration for the age-vector MDP (Algorithm 2) operates on the state space S = {0, . . . , K}N (k)

of cardinality (K + 1)N . Computing Q0 (s) at each state requires summing over all 2N failure

subsets at cost O(N ) each, giving a per-iteration cost of O N 2N (K + 1)N and space O (K + 1)N . The time-elapsed MDP (Algorithm 1), by contrast, operates on a state space of cardinality K + 1, giving time complexity O(K 2 ) and space O(K) regardless of the fleet size N . Table 1 summarizes the time and space complexity of all methods (here I is the number of iterations taken to converge). Despite not having a formal regret guarantee comparable to those of the bandit methods (see Section 4), the Kaplan-Meier approach can exhibit strong empirical performance due to its efficient use of censored data for explicit distribution learning. Table 1: Computational and space complexity of block replacement algorithms. Algorithm All Four Bandit Algorithms Kaplan-Meier Renewal Time-Elapsed MDP Age-vector MDP

Time O(KT ) O(K 2 T ) O(K 2 ) N O(N 2 (K + 1)N · I)

Space O(K) O(K) O(K) O((K + 1)N )

The choice of algorithm depends on the specific problem constraints: while bandit-based approaches offer significant computational advantages for large-scale online learning (T ≫ K), the KaplanMeier method remains competitive for moderate horizons where its superior sample efficiency outweighs its relative computational cost.

5.

Numerical Experiments

We evaluate all five data-driven algorithms (Algorithms 3 – 7) on two lifetime distributions under a common parameter regime and compare it against the age-vector MDP benchmark (Algorithm 2). Throughout this section we fix N = 2 machines, K = 12, cb = 1.0, cf = 2.6, and consider lifetimes drawn from (i) L ∼ 1 + Binomial(10, 0.5), with support {1, . . . , 11} and mean µ = 6; (ii) L ∼ 1 + Poisson(4), with support {1, 2, . . .} and mean µ = 5. Note that cb = 1 < 2cf = 5.2, so block replacement is strictly cheaper per machine than two simultaneous individual replacements, satisfying the assumption in Section 3. The binomial distribution has bounded support, while the Poisson distribution has unbounded support; for the latter, the true renewal function M (k) and cost curve c(k) are computed using a PMF truncated at Krenew = 45, which captures at least 99.9% of the probability mass. Since cf and N enter the cost function only through their product cf N , our two-machine setup 19

(which is computationally simple) with cf = 2.6 is equivalent to an N -machine fleet with per-failure cost 5.2/N , so the results scale directly to larger fleets. Our results therefore apply directly to larger fleets operating under a proportionally smaller per-failure cost. As long as cb and cf · N are fixed, the true underlying cost structure will not change, but the MDP formulation will, as its state space grows as (K + 1)N .

5.1

Age-vector MDP Formulation (Algorithm 2)

We implement Algorithm 2 for the two-machine MDP (N = 2) to illustrate the structural properties ∗ against which established in Theorem 3.4 and to compute the benchmark optimal average cost gage

the bandit and Kaplan–Meier algorithms are compared. The state space S = {0, 1, . . . , 12}2 (with |S| = 169 states), where each coordinate records the current age of a machine capped at K = 12. Value iteration is run with span-seminorm convergence tolerance ε = 10−8 . The resulting optimal policy is summarised in Table 2 via its threshold function (Theorem 3.4(d)). Table 2: Threshold function τ ∗ (b) for the optimal MDP policy under Algorithm 2. By Theorem 3.4(c), the policy is symmetric, so τ1∗ (b) = τ2∗ (a) =: τ ∗ (b): block replacement is triggered whenever either machine’s age reaches or exceeds τ ∗ (b), where b is the other machine’s current age. A threshold of 0 indicates that block replacement is triggered immediately regardless of the focal machine’s age. Other machine’s age b

0

1

2

3

4

5

6

7

8

9

10

11

12

L ∼ 1 + Binomial(10, 0.5) L ∼ 1 + Poisson(4)

4 3

4 2

3 1

2 0

0 0

0 0

0 0

0 0

0 0

0 0

0 0

0 0

0 0

∗ ∗ The optimal average costs are gage ≈ 0.3195 for Binomial and gage ≈ 0.4886 for Poisson. By

Theorem 3.4(c), the policy is symmetric in the two machine ages, so a single threshold function suffices: block replacement is triggered for machine i as soon as its age reaches τ ∗ (b), where b is the current age of the other machines. The threshold is non-increasing in b (policy monotonicity), the function is identical for both machines (symmetry), and τ ∗ (b) = 0 once the other machine is old enough (threshold structure collapses to immediate replacement). The values of τ ∗ (b) for both the Binomial and Poisson case are presented in Table 2. Note that τ ∗ = 0, meaning any age triggers block replacement. In both cases the threshold drops sharply to zero once the companion machine reaches a moderate age, and the high ratio cf N/cb = 5.2 makes proactive replacement attractive even for relatively young machines.

5.2

Bandit and Survival Analysis based methods (Algorithms 3 - 7)

We evaluate the four bandit algorithms (Algorithms 3–6) and the Kaplan–Meier renewal algorithm (Algorithm 7) under the same cost parameters and lifetime distributions as described above. Each algorithm operates over T = 10000 decision epochs, using the same sequence of underlying random failure times for all algorithms for a fair comparison. 20

The true cost curves c(k) for k ∈ [K], shown in Figure 1, reveal the following optimal block replacement intervals: (i) k ∗ = 3 for L ∼ 1 + Binomial(10, 0.5), with c(k ∗ ) ≈ 0.4282; (ii) k ∗ = 2 for L ∼ 1 + Poisson(4), with c(k ∗ ) ≈ 0.7390.

Figure 1: True long-run average cost c(k) vs k under each lifetime distribution. The optimal interval k ∗ is marked with a red diamond.

For each distribution the true optimal interval k ∗ and cost c(k ∗ ) are computed via the renewal equation (1) using a PMF truncated at Krenew = 45 to ensure accuracy under the Poisson tail. We have tested all five algorithms using simulations. The simulations need to run over long horizons for the long-run regret asymptotics to manifest itself. All five algorithms use a common random seed. The results are from a single run.

Figure 2: Cumulative regret Reg(T ) for Algorithms 3–7 over T = 10, 000 decision epochs (cb = 1, cf = 2.6, N = 2, K = 12).

Figure 2 reports the per-cycle cumulative regret (Equation 2) and Figure 3 the total-cost cumulative regret (Equation 3) over T = 10, 000 epochs across the five algorithms and two distributions. All four bandit algorithms exhibit sub-linear regret growth, consistent with the established upper bounds for Algorithms 3–6. Algorithms 5 and 6 (correlated arm bandits) accumulate substantially less regret than Algorithms 3 and 4 (independent arm bandits), confirming the theoretical O((K − k ∗ ) log T ) versus O((K − 1) log T ) separation.

Algorithms 5 and 6 exploit empirical variance

information, yielding tighter confidence bounds and faster arm elimination when the reward variance 21

g ) for Algorithms 3–7 over T = 10, 000 decision Figure 3: Cumulative regret Reg(T epochs (cb = 1, cf = 2.6, N = 2, K = 12).

differs across arms. The cumulative regret plots reveal a clear performance ordering: among the bandit algorithms, the correlated Bernstein LCB (Algorithm 6) achieves the lowest regret, consistent with its O((K − k ∗ ) log T ) guarantee and tighter variance-aware confidence bounds; the Kaplan–Meier algorithm (Algorithm 7) outperforms all bandit methods overall. By explicitly learning F from censored observations and computing c(k) via the renewal equation, the KM algorithm identifies k ∗ rapidly after the cold-start phase and thereafter incurs near-zero incremental regret. Since operating a system for T = 10,000 cycles to achieve asymptotic performance is rarely feasible in practice, we examine arm-selection behaviour of these two best-performing algorithms over the first T = 50 decision epochs; Figure 4 shows the frequency of each chosen interval k over this short horizon. The concentration of the correlated Bernstein algorithm on large arms reflects necessary exploration: the O((K −k ∗ ) log T ) regret guarantee requires every suboptimal arm above k ∗ to be pulled O(log T ) times to keep confidence bounds calibrated, and at T = 50 this amounts to visiting each such arm roughly log(50) ≈ 3.9 times. This cost is bounded and vanishes as a fraction of total cost as T → ∞. The KM algorithm closes in on the optimal arm surprisingly quickly, using the optimal arm most of the time. We find this insight most revealing and useful in practice.

5.3

Comparison of Data-Driven algorithms against MDP benchmark

The MDP formulation in Section 3.4 optimises over all age-dependent Markovian policies using complete knowledge of the lifetime distribution F , providing a lower bound on the long-run average cost. The data-driven algorithms of Section 4.1 and 4.2 operate without knowledge of F and are ∗ restricted to block replacement policies. The gap between their empirical costs and gage can be ∗ ∗ ∗ ∗ ∗ ∗ ∗ b b b decomposed as c(k (T )) − gage = c(k ) − gage + c(k (T )) − c(k ), where k (T ) := modet∈[T −W,T ] kt {z } {z } | | δlearn (T )

δstruct

is the arm most frequently selected by the algorithm in the terminal window [T − W, T ] with W = 1000 cycles, δstruct is the structural gap incurred by restricting to block policies, and δlearn (T ) is the learning gap due to not knowing F . The intermediate value c(k ∗ ) is the optimal cost of the 22

Figure 4: Histogram of replacement intervals k selected by the Correlated Bernstein LCB (Algorithm 6) and the Kaplan–Meier Renewal (Algorithm 7) algorithms over the first T = 50 decision epochs, under the Binomial(10, 0.5) (left) and Poisson(4) (right) lifetime distributions. The dashed vertical line marks the optimal interval k ∗ .

time-elapsed MDP (Section 3.3), computed analytically in O(K 2 ) operations via Algorithm 1; no simulation of the time-elapsed MDP is required. All five algorithms identify b k ∗ (T ) = k ∗ by the end of the horizon under both distributions, yielding δlearn (T ) = 0% in every case. The costs at termination are c(b k ∗ ) = c(k ∗ ) ≈ 0.4646 (Binomial) ∗ is structural. and c(b k ∗ ) = c(k ∗ ) ≈ 0.8309 (Poisson), so the entire suboptimality relative to gage ∗ ∗ ) for Binomial and δ Specifically, δstruct ≈ 0.1166 (33.5% of gage struct ≈ 0.2780 (50.3% of gage ) for

Poisson. Even an oracle with full knowledge of F that always selects the best fixed-interval block policy pays a cost above the best age-dependent policy.

6.

Conclusion

We have developed and analysed a family of data-driven algorithms for block replacement scheduling of N identical machines with unknown lifetime distribution F . The independent and correlated LCB bandit algorithms (Algorithms 3–6) achieve O(log T ) regret, matching the Lai and Robbins [1985] lower bound. The correlated variants require only O(1) direct pulls of arms k < k ∗ by exploiting the nested observation structure unique to block replacement. The Kaplan–Meier renewal algorithm (Algorithm 7) explicitly estimates F from the censored operational record and identifies k ∗ with near-zero incremental regret after a short exploration phase, empirically dominating all bandit algorithms at long horizons. At short horizons of T = 50 decision epochs, the survival analysis based algorithm already concentrates its selections near k ∗ , while the correlated Bernstein 23

algorithm spreads pulls across larger arms due to its mandatory exploration phase, a transient cost that vanishes as T → ∞. The underlying average-cost MDP admits a monotone threshold policy ∗ . under IFR lifetime distributions (Theorem 3.4), providing a gold-standard benchmark gage

The framework assumes stationary IID machine lifetimes; real systems may exhibit load-dependent degradation or temporal trends. Establishing a regret bound for Algorithm 7 and developing a Thompson sampling variant that places a prior on F and updates via censored observations are natural extensions of the proposed approaches. Bayesian approaches, including Thompson sampling (Thompson [1933], Russo et al. [2018]) and related posterior sampling methods, which places a prior on the arm reward distributions and selects arms by posterior sampling, offer competitive empirical performance in standard bandit settings and represent a natural alternative to the LCB framework developed here. Extending the framework to lifetimes that depend on environmental covariates or accumulated degradation, for example by incorporating a Cox proportional hazards model into the KM step, or by adopting contextual bandit algorithms, would substantially broaden the practical applicability of the proposed framework. ∗ The structural gap δstruct = c(k ∗ ) − gage arises because fixed-interval block replacement discards

the age information available at each decision point. A natural intermediate policy is the (k ∗ , ℓ∗ ) hybrid policy: block-replace only machines with ages at least ℓ∗ every k ∗ periods. A data-driven algorithm must learn both k ∗ and ℓ∗ from censored observations. The reduction in δstruct achievable by the optimal ℓ∗ quantifies the value of incorporating simple age-dependent information into an otherwise schedule-driven policy. The present framework assumes all N machines are statistically identical. When machines have different age profiles or lifetime distributions, as arises in mixed-age fleets following partial replacements, the bandit formulation must be extended to account for the differing contribution of each machine to the block replacement cost, a setting with closer connections to the restless bandit literature.

24

Appendix 1: Proofs of Theorems and Propositions A.

Proof of Theorem 3.4

(a) We prove by induction that the iterates {w(n) } satisfy two properties simultaneously: (A) Monotonicity: For every s ∈ S and every i ∈ [N ] with si ≤ K − 1, w(n) (s + ei ) ≥ w(n) (s). (B) Lipschitz bound: For every r ∈ S and i ∈ [N ], w(n) (r) − w(n) (r{i} ) ≤ cf , where r{i} denotes r with coordinate i reset to 0. The restriction si ≤ K − 1 in (A) is necessary to ensure s + ei ∈ S. Base case. Both properties hold trivially for w(0) ≡ 0. Inductive step. Assume (A) and (B) hold for w(n) for some n ∈ N. Fix i ∈ [N ]. For ′

brevity, write VA+′ := cf + w(n) (s[A ∪{i}] ) for the value following machine i’s failure, and ′

VA0′ := w(n) (s[A ] ) for the value following its survival. Define ΦA′ (s) := h(si +1) VA+′ + h̄(si +1) VA0′ . Conditioning on which subset A′ of the other N − 1 machines fail, (n)

Q0 (s) = cf

X

h(sj +1)

(13)

j̸=i

X

+

ωA′ (s) ΦA′ (s).

A′ ⊆[N ]\{i}

where ωA′ (s) :=

Q

j j∈A′ h(s +1)

Q

j j ∈A / ′ , j̸=i h̄(s +1) ≥ 0 with ′

P

′

ωA′ (s) = 1. In s[A ] , machine

A′

i survives and its age advances to si + 1; in s[A ∪{i}] , machine i fails and resets to age 0. Proof of (A) for w(n+1) . Fix s ∈ S with si ≤ K − 1, so that s + ei ∈ S. Replace the lead-in and equation (14) with: ′

Let se[A ] denote the next state from s + ei when machines in A′ fail and machine i survives: ′

se[A ] has coordinate i equal to si + 2 (with zero weight when si + 1 = K), coordinate j = 0 ′

′

for j ∈ A′ , and coordinate j = sj + 1 for j ∈ / A′ ∪ {i}. Note se[A ] ⪰ s[A ] coordinatewise and ′

′

(e s[A ] ){i} = s[A ∪{i}] . Note that δhi := h(si +2) − h(si +1) ≥ 0 by the IFR assumption. Define ′

′

′

′

R1 (A′ ) := cf + w(n) (s[A ∪{i}] ) − w(n) (e s[A ] ) and R2 (A′ ) := w(n) (e s[A ] ) − w(n) (s[A ] ). Then (n)

(n)

Q0 (s + ei ) − Q0 (s) X   = ωA′ (s) δhi · R1 (A′ ) + h̄(si +1) · R2 (A′ ) . A′

25

(14)

′

′

The term R2 (A′ ) ≥ 0 by inductive hypothesis (A), since se[A ] ⪰ s[A ] . For R1 (A′ ), applying (B) ′

′

′

′

′

to se[A ] at coordinate i gives w(n) (e s[A ] ) − w(n) ((e s[A ] ){i} ) = w(n) (e s[A ] ) − w(n) (s[A ∪{i}] ) ≤ cf , (n)

(n)

so R1 (A′ ) ≥ 0. Since δhi ≥ 0, both terms are non-negative, giving Q0 (s + ei ) ≥ Q0 (s). (n)

(n)

Proof of (B) for w(n+1) . Fix r ∈ S and i ∈ [N ]. Since Q0 (r) ≥ Q0 (r{i} ) by (A) (as r{i} + ei has coordinate i equal to 1 ≤ K), and min(a, c) − min(b, c) ≤ a − b for any reals a ≥ b: (n)

(n)

w(n+1) (r) − w(n+1) (r{i} ) ≤ Q0 (r) − Q0 (r{i} ). We apply (13) to both r and r{i} and subtract. When machine i fails, both r and r{i} produce ′

′

next-state coordinate i = 0, so r[A ∪{i}] = (r{i} )[A ∪{i}] . When machine i survives from r it ′

′

ages to ri + 1; from r{i} it ages to 1. Define B1 (A′ ) := cf + w(n) (r[A ∪{i}] ) − w(n) ((r{i} )[A ] ) ′

′

and B2 (A′ ) := w(n) (r[A ] ) − w(n) ((r{i} )[A ] ). ′

′

′

′

By inductive hypothesis (A), r[A ] ⪰ (r{i} )[A ] at coordinate i and (r{i} )[A ] ⪰ r[A ∪{i}] at ′

coordinate i, so both B1 (A′ ) ≥ 0 and B2 (A′ ) ≥ 0. Applying (B) to (r{i} )[A ] at coordinate ′

i gives B1 (A′ ) ≤ cf , and applying (B) to r[A ] at coordinate i gives B2 (A′ ) ≤ cf . Since h(ri + 1) − h(1) ≥ 0 by IFR, (n)

(n)

Q0 (r) − Q0 (r{i} ) X 
= ωA′ (r) h(ri +1) − h(1) B1 (A′ ) A′

 +h̄(ri +1) B2 (A′ ) . Using

P

A′ ωA′ = 1, h̄(r

(15)

i + 1) ≤ 1, and the bounds B , B ≤ c , we obtain 1 2 f

(n) (n) Q0 (r) − Q0 (r{i} ) ≤ cf h(ri +1) − h(1) + 1 − h(ri +1) = cf (1 − h(1)) ≤ cf , establishing (B) for n + 1. Transfer to w∗ . Under the unichain condition (Theorem 3.3), sp(w(n+1) − w(n) ) → 0. Since w(n) (0) = 0 for all n, the sequence {w(n) } is uniformly bounded and converges pointwise to w∗ . Properties (A) and (B) are preserved under pointwise limits, so w∗ satisfies: (A) w∗ (s + ei ) ≥ w∗ (s) for all s with si ≤ K − 1, i.e., w∗ is non-decreasing in each coordinate. (B) w∗ (r) − w∗ (r{i} ) ≤ cf for all r ∈ S, i ∈ [N ]. This completes the proof of part (a). (b) Define the advantage function ∆(s) := Q∗0 (s)−Q∗1 (s), where Q∗1 (s) = cb +h∗ (0) is constant in s. The optimal policy satisfies π ∗ (s) = 1 if and only if ∆(s) ≥ 0. Since Q∗0 (s + ei ) ≥ Q∗0 (s) for all i ∈ [N ], we have Q∗0 (s) is non-decreasing in each coordinate, and hence ∆(s) is non-decreasing in each coordinate. If π ∗ (s) = 1 then ∆(s) ≥ 0, and for any s′ ≥ s, ∆(s′ ) ≥ ∆(s) ≥ 0, so 26

π ∗ (s′ ) = 1. (c) By symmetry of the problem setup, permuting machines 1, 2, · · · , N yields an equivalent MDP. Hence both the value function and optimal policy inherit this symmetry. (d) Fix i ∈ [N ] and s−i ∈ {0, 1, · · · , K}N −1 . Define τi∗ (s−i ) := min{si ∈ {0, . . . , K} : π ∗ (si , s−i ) = 1}, setting τi∗ (s−i ) = ∞ if no such si exists. If si ≥ τi∗ (s−i ), then (si , s−i ) ≥ (τi∗ (s−i ), s−i ) coordinatewise and π ∗ (τi∗ (s−i ), s−i ) = 1 by definition, so Theorem 3.4(b) gives π ∗ (si , s−i ) = 1. Conversely, if si < τi∗ (s−i ), then π ∗ (si , s−i ) = 0 by minimality of τi∗ (s−i ). This establishes the threshold representation. For monotonicity, suppose s′−i ≥ s−i componentwise and τi∗ (s−i ) < ∞. Setting si0 := τi∗ (s−i ), we have π ∗ (si0 , s−i ) = 1 and (si0 , s′−i ) ≥ (si0 , s−i ) coordinatewise, so Theorem 3.4(b) gives π ∗ (si0 , s′−i ) = 1, hence τi∗ (s′−i ) ≤ si0 = τi∗ (s−i ).

B.

■

Proof of Theorem 4.3

Using the regret decomposition in Equation 2, it suffices to show that E[nk (T )] = O(1) for all k < k ∗ . Fix k < k ∗ . Since k ∗ is assumed to be unique, ∆k = c(k) − c(k ∗ ) > 0. Define the critical 2 sample size n0 = ⌈16b2 ∆−2 k log(KT )⌉ where b = cf N . By Hoeffding’s inequality, for any j ∈ [K]

with mj (t) ≥ n0 and mk∗ (t) ≥ n0 ,  P Define the event Gt =

(m)

b cj

− c(j) ≥ ∆4k



  m∆2 ≤ 2 exp − 8b2k ≤

2 . KT 2

n o n o (m) (m) b ck (t) ≤ c(k) + ∆k /4 ∩ b ck∗ (t) ≥ c(k ∗ ) − ∆k /4 . By union bound,

P(Gtc ) ≤ 4/(KT 2 ). For m ≥ n0 , r q q 2) 2 log(KT 2 ) 2 log t = b 16b22∆log(KT = ∆4k . b −2 m ≤b n0 log(KT 2 ) k

∆k ∆k ∆k On Gt , the LCB values therefore satisfy (i) LCBCH k (t) ≥ c(k) − 4 − 4 = c(k) − 2 and also (ii) ∆k ∆k ∆k ∗ ∗ ∗ LCBCH k∗ (t) ≤ c(k ) + 4 + 4 = c(k ) + 2 . Since c(k) − c(k ) = ∆k , by (i) and (ii), we have ∆k CH ∗ LCBCH k (t) ≥ c(k ) + 2 > LCBk∗ (t),

so arm k is not selected at epoch t. Arm k can only be pulled when mk (t) < n0 , or mk∗ (t) < n0 , or Gtc occurs. Thus, E[nk (T )] ≤ 2n0 +

T X

∞

P(Gtc ) ≤ 2n0 +

t=1

4 X 1 2π 2 = 2n0 + = O(1). 2 K t 3K t=1

27

C.

Proof of Theorem 4.4

As in the proof of the previous theorem, it suffices to show that E[n k < k ∗. k (T )] = O(1) for all 

Fix k < k ∗ . Since k ∗ is unique, ∆k = c(k) − c(k ∗ ) > 0. Let σk2 = Var

C (k) k

and σk2∗ = Var

∗

C (k ) k∗

denote the variances, with σk2 , σk2∗ ≤ b2 where b = cf N . om l n 2 om l n 2 96σk∗ 96σ 2) 2 ), 48b log(KT 2 ) ∗ log(KT and n = max Define n0 = max ∆2k log(KT 2 ), 48b log(KT 0,k ∆k ∆k ∆2 k

k

to be the critical sample sizes. By Bernstein’s inequality, for m ≥ n0 and any arm j,  P

∆k (m) b cj − c(j) ≥ 4



m∆2 /16 ≤ 2 exp − 2 k 2σj + b∆k /6

! ≤

2 . KT 2

  (m) 2 Moreover, for m ≥ n0 , by Bernstein’s inequality, P (b σj )2 ≥ 2σj2 ≤ KT 2. Define n o n o o n o n (m) (m) (m) (m) Gt := |b ck (t)−c(k)| ≤ ∆4k ∩ |b σk∗ )2 (t) ≤ 2σk2∗ . ck∗ (t)−c(k ∗ )| ≤ ∆4k ∩ (b σk )2 (t) ≤ 2σk2 ∩ (b (m)

When mk (t) ≥ n0 and mk∗ (t) ≥ n0,k∗ , the union bound gives P(Gtc ) ≤ 8/(KT 2 ). Since (b σk∗ )2 (t) ≤ 2σk2∗ on Gt and mk∗ (t) ≥ n0,k∗ ≥ 2, we have ∆k ∗ LCBCB − k∗ (t) ≤ c(k ) + 4 By definition of n0,k∗ ,

q

4.8σk2∗ log t/n0,k∗ ≤

s

2.4 · 2σk2∗ log t 3.6b log t − . n0,k∗ n0,k∗

q √ 4.8σk2∗ /96σk2∗ · ∆k = ∆k / 20, and 3.6b log t/n0,k∗ ≤

√ 3.6b/(48b) · ∆k = ∆k / 50 (both uniformly in t ≤ T ). Hence ∗ LCBCB k∗ (t) ≤ c(k ) +

∆k ∆k ∆k ∆k ∆k ∆k − √ − √ ≤ c(k ∗ ) + − < c(k ∗ ) + . 4 4 5 8 20 50

  (m) 2 (m) Since σ bk (t) ≤ 2σk2 ≤ b2 on Gt and b ck (t) ≥ c(k) − ∆k /4, we have LCBCB k (t) ≥ c(k) −

∆k − 4

s

2.4 · 2σk2 log t 3.6b log t − . mk (t) mk (t)

√ Even if mk (t) = 1 (Hoeffding bound applies), LCBCB k (t) ≥ c(k) − ∆k /4 − b 2 log t, and arm k is not selected once n0 is reached regardless. For mk (t) ≥ n0 , a similar argument gives LCBCB k (t) ≥ c(k) −

∆k ∆k ∆k ∆k − > c(k) − = c(k ∗ ) + . 4 5 2 2

CB ∗ ∗ Thus, we get LCBCB k∗ (t) < c(k ) + ∆k /8 < c(k ) + ∆k /2 < LCBk (t), so arm k is not selected.

28

Thus arm k can only be pulled when mk (t) < n0 , or mk∗ (t) < n0,k∗ , or Gtc holds. Hence E[nk (T )] ≤ n0 + n0,k∗ +

T X

P(Gtc ) ≤ n0 + n0,k∗ +

t=1

8 π2 · = O(1). K 6 ■

D.

Proof of Theorem 4.5

Conditional on the chosen policy kt and the failure history up to τ j,tj , the residual lifetime Lj,tj Yt

Yt +1

has

the same distribution F as the original lifetime, and is independent of the deterministic quantity C j,t = kt − τ j,tj . This corresponds to Type I censoring with a fixed censoring time, which trivially Yt

satisfies the independent censoring assumption. Thus the uniform strong consistency follows from the Glivenko-Cantelli theorem for randomly right-censored data (see Breslow and Crowley [1974] a.e. and Gill [1983]), and hence supt∈[K] |SbT (t) − S(t)| −−→ 0. Since FbT (t) = 1 − SbT (t), uniform strong consistency of FbT follows immediately.

E.

■

Proof of Theorem 4.6 a.s.

cT , We first show that b cT (k) −−→ c(k) as T → ∞ for every k ∈ [K]. For pointwise convergence of M cT (0) = 0 = M (0) holds by definition. Suppose we proceed by induction on k. The base case M a.s. cT (j) − M −→ M (j) for all j < k. The estimated renewal function satisfies

cT (k) = FbT (k) + M

k X

cT (k − j). fbT (j) M

j=1 a.s. By Theorem 4.5, FbT (k) −−→ F (k), so a.s. fbT (j) = FbT (j) − FbT (j − 1) −−→ f (j). a.s. cT (k − j) − Combined with the inductive hypothesis M −→ M (k − j), each of the finitely many terms a.s. c converges almost surely, giving MT (k) −−→ M (k). Thus for each k ∈ [K],

b cT (k) =

cT (k) a.s. cb + cf N M (k) cb + cf N M −−→ = c(k). k k

Since k ∗ is assumed to be unique, δ := mink̸=k∗ ∆k > 0. For all sufficiently large nT , |b cT (k)−c(k)| < ∗ ∗ ∗ ∗ δ/2 for every k ∈ [K], whence b cT (k) > b cT (k ) for all k ̸= k and thus b k = k . The almost-sure T

cost convergence |b cT (b kT∗ ) − c(k ∗ )| → 0 follows immediately. 29

■

Appendix 2: Algorithm Pseudocodes This appendix collects the pseudocodes for all seven algorithms developed in the main text. Algorithms 1 and 2 are the value iteration procedures for the two MDP formulations of Section 3 and require knowledge of the lifetime distribution F ; they serve as benchmarks rather than datadriven policies. Algorithms 3 - 6 are the bandit-based online learning algorithms of Section 4.1, which operate without knowledge of F and differ in whether arms are treated as independent or correlated and whether Hoeffding or Bernstein concentration inequalities are used to construct confidence bounds. Algorithm 7 is the Kaplan–Meier renewal algorithm of Section 4.2, which explicitly estimates F from censored operational data. Table 1 of the main text summarises the time and space complexity of all seven algorithms.

A2.1

Value Iteration for the Time-Elapsed MDP

Algorithm 1 implements value iteration for the time-elapsed MDP of Section 3.3. The state space has cardinality K + 1, independent of fleet size N , making this procedure computationally inexpensive. The input requires the lifetime PMF f , from which the renewal function M (·) and renewal mass function m(·) are computed via the discrete renewal equation (1). Convergence is monitored via the span seminorm sp(d(n) ) := maxk d(n) (k) − mink d(n) (k) of successive iterates, which contracts to zero under the unichain condition established in Theorem 3.1. Upon termination, the midpoint estimate b l∗ is ε/2-accurate and the extracted policy π b∗ is ε/2-optimal (Theorem 3.1).

The output threshold ks = k ∗ established in Theorem 3.2 directly grounds the bandit formulation of Section 4.1. Algorithm 1: Value Iteration for the Time-Elapsed MDP Input: PMF f , parameters K, N, cf , cb , tolerance ε > 0 Output: b l∗ , w(n) , π b∗ Compute M (k) for k = 1, . . . , K via the renewal equation (1); Set m(k) ← M (k) − M (k − 1) for k = 1, . . . , K; Initialize w(0) (k) ← 0 for all k ∈ {0, 1, . . . , K}; repeat for k = 1, . . . , K − 1 do
u(k) ← min cf N m(k) + w(n) (k + 1), cb + w(n) (0) ; u(K) ← cb + w(n) (0); w(n+1) (k) ← u(k) − u(0) for all k ; // Normalize at reference state (n) (n+1) (n) d (k) ← w (k) − w (k) for all k; until sp(d(n) ) < ε;
b l∗ ← u(0) + 12 mink d(n) (k) + maxk d(n) (k) ; π b∗ (k) ← arg mina∈{0,1} Qa (k) for all k, where Q0 (k) = cf N m(k) + w(n) (k + 1), Q1 (k) = cb + w(n) (0);

30

A2.2

Value Iteration for the Age-Vector MDP

Algorithm 2 implements value iteration for the age-vector MDP of Section 3.4. Unlike Algorithm 1, the state space S = {0, 1, . . . , K}N has cardinality (K + 1)N , so the procedure is tractable only for small fleet sizes N ; in the numerical experiments of Section 5 we use N = 2, giving |S| = 169 states. The per-iteration cost is O(N 2N (K + 1)N ) because computing Q0 (s) at each state requires summing over all 2N failure subsets. Convergence is again monitored via the span seminorm, and the termination guarantees of Theorem 3.3 apply. The output optimal average cost gb∗ and policy π b∗

∗ against which all data-driven algorithms are evaluated serve as the gold-standard benchmark gage

in Section 5.3. By the symmetry established in Theorem 3.4, the optimal policy depends only on
threshold values. the multiset of machine ages, reducing effective storage from (K + 1)N to K+N N

Algorithm 2: Value Iteration for Average-Cost MDP Input: CDF F , costs (cb , cf ), N , K, tolerance ε > 0 Output: Optimal average cost gb, value function w(k) , policy π b∗ w(0) ≡ 0; k ← 0; sref ← 0 repeat foreach s ∈ S do  P P Q Q j +1) j +1) w (k) (s[A] ) Q0 (s) ← cf h(sj +1) + h(s h̄(s j∈A j ∈A / j

A⊆[N ]

v (k+1) (s) ← min{Q0 (s), cb } foreach s ∈ S do w(k+1) (s) ← v (k+1) (s) − v (k+1) (sref ) d(k) ← w(k+1) − w(k) ; sp(d(k) ) ← maxs d(k) (s) − mins d(k) (s); k ←k+1 until sp(d(k) ) < ε
gb ← v (k) (sref ) + 12 mins d(k) (s) + maxs d(k) (s) foreach s ∈ S do π b∗ (s) ← 1[cb ≤ Q0 (s)] return gb; w(k) ; π b∗

A2.3

// recompute Q0 using w(k)

Independent Hoeffding LCB

Algorithm 3 is the simplest of the four bandit algorithms. Each arm k ∈ [K] is treated as an independent learning problem, and the lower confidence bound is constructed from Hoeffding’s p inequality with exploration bonus b 2 log(t)/nk (t), where b = cf N is the half-range of the perunit-time cost. The initialization phase pulls each arm once in order of increasing k before switching to the LCB rule; this ensures every arm has at least one sample before confidence bounds are evaluated. The algorithm achieves O((K − 1) log T ) cumulative regret (Theorem 4.1), matching the Lai–Robbins lower bound up to constants. 31

Algorithm 3: Independent Hoeffding LCB Input: K, T , cb , cf , N Output: k1 , . . . , kT foreach k ∈ [K] do nk ← 0; Sk ← 0 for t = 1 to T do if ∃ k : nk = 0 then kt ← min{k : nk = 0} n o p else kt ← arg min Sk /nk − b 2 log(t)/nk k∈[K]

Observe total failures Yt over interval kt nkt ← nkt + 1; Skt ← Skt + (cb + cf Yt )/kt

A2.4

Independent Bernstein LCB

Algorithm 4 replaces the Hoeffding exploration bonus with the empirical Bernstein bound of Audibert et al. [2009], which substitutes the worst-case range b2 with the sample variance σ bk2 (t). When nk (t) = 1 the Hoeffding bonus is used as a fallback since the sample variance is undefined; for nk (t) ≥ 2 the Bernstein LCB of equation (11) applies. The algorithm achieves the same O((K − 1) log T ) regret rate as Algorithm 3 (Theorem 4.2) but with tighter instance-dependent constants when arm reward variances differ. Algorithm 4: Independent Bernstein LCB Input: K, T , cb , cf , N Output: k1 , . . . , kT foreach k ∈ [K] do nk ← 0; Sk ← 0; Qk ← 0 for t = 1 to T do if ∃ k : nk = 0 then kt ← min{k : nk = 0} else foreach k ∈ [K] do b ck ← Sk /nk ; σ bk2 ← (Qk − Sk2 /nk )/(nk − 1) if nk ≥ 2, else b2 kt ← arg min LCBIB // Eq. (11); Hoeffding bound used when nk = 1 k (t) k∈[K]

Observe Yt ;

A2.5

ct ← (cb + cf Yt )/kt ;

nkt ← nkt + 1;

Skt ← Skt + ct ;

Qkt ← Qkt + c2t

Correlated Hoeffding LCB

Algorithm 5 exploits the nested observation structure of Proposition 4.1: pulling arm k can be used to infer the costs of all shorter intervals j ≤ k can be reconstructed at no additional P cost. Consequently, the effective sample count for arm j is mj (t) = K l=j nl (t) ≥ nj (t), and the LCB is computed using mj (t) in place of nj (t). The algorithm achieves O((K − k ∗ ) log T ) regret (Theorem 4.3), a strict improvement over the O((K − 1) log T ) bound of Algorithm 3 whenever k ∗ > 1, because arms k < k ∗ require only O(1) direct pulls. 32

Algorithm 5: Correlated Hoeffding LCB Input: K, T , cb , cf , N Output: k1 , . . . , kT foreach k ∈ [K] do mk ← 0 Sk ← 0 for t = 1 to T do if ∃ k : mk = 0 then kt ← max{k : mk = 0} n o p else kt ← arg min Sk /mk − b 2 log(t)/mk k∈[K]

Observe cumulative failures: Zj ←

j P

Yi , for j = 1, . . . , kt

i=1

foreach j = 1, . . . , kt do mj ← mj + 1

A2.6

Sj ← Sj + (cb + cf Zj )/j

Correlated Bernstein LCB

Algorithm 6 combines the nested observation structure of Algorithm 5 with the empirical Bernstein (m)2

confidence bound of Algorithm 4. The sample variance σ bk

(t) is now computed from the mk (t)

correlated samples rather than the nk (t) direct pulls. The algorithm achieves O((K − k ∗ ) log T ) regret (Theorem 4.4), with the same improvement over independent arms as Algorithm 5 but potentially tighter constants due to variance adaptation.

Algorithm 6: Correlated Bernstein LCB Input: K, T , cb , cf , N Output: k1 , . . . , kT foreach k ∈ [K] do mk ← 0 Sk ← 0 Qk ← 0 for t = 1 to T do if ∃ k : mk = 0 then kt ← max{k : mk = 0} else foreach k ∈ [K] do σ bk2 ← (Qk − Sk2 /mk )/(mk − 1) if mk ≥ 2, else b2 // From Eq. (12) kt ← arg min LCBCB k (t) k∈[K]

Set C (j) ← (cb + cf

j P

Yi )/j for j = 1, . . . , kt

i=1

foreach j = 1, . . . , kt do mj ← mj + 1

A2.7

Sj ← Sj + C (j)

Qj ← Qj + (C (j) )2

Kaplan–Meier Renewal Algorithm

Algorithm 7 takes a fundamentally different approach from the bandit algorithms: rather than directly estimating the per-unit-time cost c(k) for each arm, it estimates the lifetime distribution F nonparametrically from the accumulated censored observations and uses the estimated Fb to compute b c(k) via the renewal equation 1. At each block replacement, every machine that has not failed contributes a right-censored observation, and every complete lifetime between failures b of Theorem 4.5 correctly contributes an uncensored observation; the Kaplan–Meier estimator S(t) incorporates both types. The algorithm has two hyperparameters: the exploration rate ε governing 33

√ the decaying exploration schedule εt = ε/ t, and the refit frequency R, which controls how often the KM estimator and renewal function are recomputed from the accumulated dataset DT . A cold-start phase (t ≤ 5) ensures that all K arms are pulled at least once before exploitation begins, providing the KM estimator with observations before the first policy update. Despite the absence of a formal regret guarantee, the algorithm identifies k ∗ rapidly after the cold-start phase and thereafter incurs near-zero incremental regret, as established by Theorem 4.6 and confirmed in Section 5.2. Algorithm 7: Kaplan–Meier Renewal Algorithm Input: (cb , cf ), N , Tmax , refit frequency R, ε ∈ (0, 1), K Output: b k ∗ , {b c(k)}K k=1 b D ← ∅; t ← 0; k ∗ ← Uniform([K]) while t < Tmax do t←t+1 √ if t ≤ 5 or Bernoulli(ε/ t) = 1 then kt ← Uniform([K]) else kt ← b k∗ for j = 1 to N do τ ←0 while τ < kt do sample ℓ ∼ F if τ + ℓ ≤ kt then D ← D ∪ {(ℓ, 1)}; τ ←τ +ℓ else D ← D ∪ {(kt − τ, 0)}; break if t mod R = 1 or t = 6 then Compute Fb via KM on D;

b − 1) − S(s); b c(0) ← 0 fb(s) ← S(s M s P c(s) ← c(s − r)] foreach s = 1, . . . , K do M fb(r)[1 + M r=1

c(k))/k foreach k = 1, . . . , K do b c(k) ← (cb + cf N M b k ∗ ← arg min b c(k) k∈[K]

return b k∗ ;

{b c(k)}K k=1

34

References Andersen JF, Andersen AR, Kulahci M, Nielsen BF (2022) A numerical study of Markov decision process algorithms for multi-component replacement problems. European Journal of Operational Research 299:898–909. Audibert J-Y, Munos R, Szepesvári C (2009) Exploration–exploitation tradeoff using variance estimates in multi-armed bandits. Theoretical Computer Science 410(19):1876–1902. Auer P, Cesa-Bianchi N, Fischer P (2002) Finite-time analysis of the multiarmed bandit problem. Machine Learning 47:235–256. Barlow RE, Hunter L (1960) Optimum preventive maintenance policies. Operations Research 8(1):90–100. Baxter LA, Li L (1995) Nonparametric confidence intervals for the renewal function with censored data. Journal of Nonparametric Statistics 4:317–326. Berg M (1976) Optimal replacement policies for two-unit machines with increasing running costs. Stochastic Processes and their Applications 4:89–106. Breslow N, Crowley J (1974) A large sample study of the life table and product limit estimates under random censorship. The Annals of Statistics 2(3):437–453. Bubeck S, Cesa-Bianchi N (2012) Regret analysis of stochastic and nonstochastic multi-armed bandit problems. Foundations and Trends in Machine Learning 5:1–122. Caron S, Kveton B, Lelarge M, Bhagat S (2012) Leveraging side observations in stochastic bandits. Proceedings of the 28th Conference on Uncertainty in Artificial Intelligence 142–151. Cho DI, Parlar M (1991) A survey of maintenance models for multi-unit systems. European Journal of Operational Research 51(1):1–23. Corman F, Kraijema S, Godjevac M, Lodewijks G (2017) Optimizing preventive maintenance policy: a data-driven application for a light rail braking system. Journal of Risk and Reliability 231:534– 545. Dekker R, Wildeman RE, van der Duyn Schouten F (1997) A review of multi-component maintenance models with economic dependence. Mathematical Methods of Operations Research 45:411–435. Deprez L, Antonio K, Boute R (2021) Pricing service maintenance contracts using predictive analytics. European Journal of Operational Research 290:530–545. Ge P, Kulkarni VG (2024) Data driven machine maintenance: a renewal reward bandit framework. SSRN Electronic Journal doi:10.2139/ssrn.4959529. 35

Gill RD (1983) Large sample behaviour of the product-limit estimator on the whole line. The Annals of Statistics 11(1):49–58. Grübel R, Pitts SM (1993) Nonparametric estimation in renewal theory I: the empirical renewal function. The Annals of Statistics 21:1431–1451. Kaplan EL, Meier P (1958) Nonparametric estimation from incomplete observations. Journal of the American Statistical Association 53(282):457–481. Kulkarni VG (2009) Modeling and Analysis of Stochastic Systems, 2nd ed. (Chapman and Hall/CRC, Boca Raton). Kurt M, Kharoufeh JP (2010) Monotone optimal replacement policies for a Markovian deteriorating system in a controllable environment. Operations Research Letters 38:273–279. Kveton B, Wen Z, Ashkan A, Szepesvári C (2015) Tight regret bounds for stochastic combinatorial semi-bandits. Proceedings of the Eighteenth International Conference on Artificial Intelligence and Statistics 535–543. Lai TL, Robbins H (1985) Asymptotically efficient adaptive allocation rules. Advances in Applied Mathematics 6(1):4–22. Lattimore T, Szepesvári C (2020) Bandit Algorithms (Cambridge University Press, Cambridge). Mannor S, Shamir O (2011) From bandits to experts: on the value of side-observations. Advances in Neural Information Processing Systems 24. Puterman ML (1994) Markov Decision Processes: Discrete Stochastic Dynamic Programming (John Wiley & Sons, New York). Ross SM (1996) Stochastic Processes, 2nd ed. (John Wiley & Sons, New York). Russo DJ, Van Roy B, Kazerouni A, Osband I, Wen Z (2018) A tutorial on Thompson sampling. Foundations and Trends in Machine Learning 11:1–96. Shaked M, Zhu H (1992) Some results on block replacement policies and renewal theory. Journal of Applied Probability 29(4):932–946. Sheu S-H, Jhang J-P (1996) A generalized group maintenance policy. European Journal of Operational Research 96:232–247. Thompson WR (1933) On the likelihood that one unknown probability exceeds another in view of the evidence of two samples. Biometrika 25:285–294. Vardi Y (1982) Nonparametric estimation in renewal processes. The Annals of Statistics 10:772– 785.

36