One possible alternative approach is incorporating deep reinforcement learning (DRL) techniques. These techniques are sequential since they incorporate lookahead, and they are designed to attempt to reach an optimal solution. Moreover, a Bayesian viewpoint can be adopted through incorporating a belief state over the unknown parameters. This belief state transforms the original POMDP into a belief MDP (see second paragraph of page 2 of this paper) where the observations are now belief states (see the wiki for more formal details). The belief state is generally more informative than a simple observation of a POMDP, which should aid in training speed and stability.
These DRL methods are appealing for the posed problem setting since (i) many have off-the-shelf implementations for fast use, (ii) they natively navigate the exploration/exploitation dilemma, and (iii) they can be quickly applied to any new variant of the posed problem as long as it is modeled correctly (e.g. as a belief MDP). Since the question is asking about implementation details, this answer will detail the some aspects of the approach discussed in aforementioned paper. The first major consideration is the construction of the state for use in the neural network. Most practical DRL algorithms will optimize over a fixed planning horizon which we will assume is of length $T$. The first order of business is to create a one-hot encoding vector that indicates the current experiment number. This initial vector will be $0$ everywhere except in index $k-1$ during experiment number $k$ (starting from index $0$): $$[\underbrace{0, 0, 0, \ldots, 0, \underbrace{1}_{\mbox{index}\ k-1}, 0, \ldots, 0}_{T}].$$ This vector needs to be augmented with the designs of each experiment as follows: $$[\underbrace{x_{1, t_0}, x_{2, t_0}, x_{1, t_1}, x_{2, t_1},\ldots, x_{1, t_{k-2}}, x_{2, t_{k-2}}, x_{1, t_{k-1}}, x_{2, t_{k-1}}, 0, 0, \ldots, 0, 0}_{2T}].$$ All designs yet to be chosen can simply be set to 0. Finally, the history of observations needs to be augmented as well, with future observations set to 0: $$[\underbrace{y_0, y_1, y_2, \ldots, y_{k-2}, y_{k-1}, 0, \ldots, 0}_{T}].$$ This entire vector is a vector of length $T+2T+T = 4T$, which is definitely a reasonable input size for a neural network under a planning horizon of ~30 days. Note that the designs and observations are all that is necessary (along with the original prior distribution) to construct the belief state. Since the state vector includes these quantities as a subset, also known as the information set $$I_k = \{x_{1, t_0}, x_{2, t_0}, y_0, x_{1, t_1}, x_{2, t_1}, y_1, \ldots, x_{1, t_{k-1}}, x_{2, t_{k-1}}, y_{k-1}\},$$ then the belief state is indeed encoded in nonparametric form in the state vector. By not constructing the belief state explicitly, intermediate belief states do not require any additional computation, which generally is non-negligible (see page 10 near equation 26 in the aforementioned paper for other details and discussion of a nonparametric belief state).
Following the state construction of the belief MDP, the issue of action constraints needs to be addressed since there is a design budget in the posed problem. Some off-the-shelf DRL libraries allow for parametric or constrained action spaces that are dependent on the state. In this case, augmenting the state vector with the remaining budget may aid in training speed due to a more explicit state representation (as opposed to the neural network needing to learn budget remaining from the designs).
After the above primary considerations are decided, all that is left is to code the belief MDP as a RL environment and apply some DRL algorithms. Some of the references in the linked question can help provide insight regarding the application of DRL to such environments. The posed problem seems non-trivial, open-ended, and may require substantial experimentation; therefore, this answer will be inherently limited in scope. Here are some more major considerations that may aid during the experimental process:
Reward: Some applications provide a reward solely for the information gain in the unknown parameters. If information gain should be incorporated into the reward function, check out the references in the link question for ideas.
Recomputing decisions: Since the entire analysis can be recomputed each night, the algorithm can be rerun with a shortened state vector. The new belief state can be computed and used as the prior distribution in the new algorithm. A shorter time horizon is generally easier for DRL algorithms and could either confirm the future design choices of past algorithms or even improve upon them.
Belief state: The belief state representation was nonparametric and somewhat primitive. A better belief state representation may aid in training efficiency (e.g. parameters of the posterior distribution, if they are available).
POMDP solutions: Since a belief MDP is a special case of POMDPs, techniques that solve POMDPs might be good baselines. For instance, this paper shows how memory-based architectures provide strong baselines for solving POMDPs with DRL.
Good luck, and keep us posted with any success stories!
Addendum: Other Belief State Representations
Since you were curious about other possible belief state representations, I will detail two general possibilities below.
The first possibility is to make use of a conjugate prior distribution (wiki). Given a likelihood function, a conjugate prior can sometimes be chosen such that the posterior distribution is from the same distribution family as the prior. This approach prevents costly computational Bayesian updates. For reference, this paper makes use of a conjugate prior in section 4.1.
The second possibility is to use a grid discretization over the unknown parameters. In the case of two unknown parameters, a grid of rectangles can be placed over the possible (or high-likelihood) locations of the unknown parameters. Each rectangle can be associated with a single set of parameters, often the parameters located at the center of the rectangle. In this manner, a Bayesian update can be performed by computing the product of the likelihood and prior at each rectangle and then normalizing by dividing by the sum of such products across all rectangles. This method is appealing due to the ability of viewing each posterior as a heatmap. The previously linked papers by Shen & Huan and Huan & Marzouk serve as good references for this method.