I am interested in potentially using PyStan for an optimization problem I am dealing with.

I have a large amount of historical sales data for six different categories of products. Each product has a size attribute. I want to optimize the size of each product, subject to the constraints that product 1 must be smaller than product 2, which must be smaller than product 3, etc.

I could use a Gaussian process for estimating the relationship between sales & size within each product category (or even just a simple regression), and then find the maximum of the resulting equation. However that wouldn’t leave my results subject to the constraints.

Is it possible to model the relationship between product category, size, and sales, subject to the constraints, and get the maximum for each in PyStan? Or am I barking up the wrong tree?