Solving High-Dimensional Multiobjective Optimization Problems

Solving high-dimensional blackbox optimization problems is very hard. Solving them in the multiobjective sense is even harder. Doing this efficiently on a limited budget could be considered a grand challenge type problem.

To quote our FAQ:

  • The key issue is that global optimization is expensive. At a fundamental level, no blackbox solver can guarantee global convergence without densely sampling the design space, which is exponentially expensive when n (number of design variables) is large. So what can you do? You can switch to using local modeling methods, whose costs generally only grow linearly in the dimension. You will not get any global convergence guarantees, but in many cases, you will still be able to approximately solve your problem.

  • If you have a lot of design variables, then you might do better with a local solver, by switching your surrogate to the LocalGaussRBF surrogate. If you are using the LBFGSB optimizer, then you will also need to switch to the TR_LBFGSB optimizer.

  • The majority of ParMOO’s overhead comes from fitting the surrogate models and solving the scalarized surrogate problems. If you followed the quickstart, then the default method for surrogate modeling was to fit a Gaussian process. For numerical stability reasons, we fit our Gaussian processes via a symmetric-eigensolve, which is not cheap. Then you may have to evaluate the Gaussian process thousands of times while solving the surrogate problem. All of this expense adds up, especially if you are using a large total budget, since the cost of fitting Gaussian processes grows cubically with the number of data points. One solution is to switch to using a LocalGaussRBF surrogate, which does not use the entire database when fitting surrogates, and therefore is more scalable for handling large budgets.

We will attempt to solve a convex 50-design variable, 2-objective problem on a budget of just 1000 simulation evaluations. Going off a modification to the quickstart, this will produce the following script. We also have a similar example in the solver_farm.


import numpy as np
import pandas as pd
from parmoo import MOOP
from parmoo.searches import LatinHypercube
from parmoo.surrogates import LocalGaussRBF
from parmoo.acquisitions import RandomConstraint, FixedWeights
from parmoo.optimizers import TR_LBFGSB
from parmoo.objectives import single_sim_out
import logging

# Fix the random seed for reproducibility
np.random.seed(0)

# Switch to using the TR_LBFGSB solver to solve surrogate problems in
# a trust region with multi-start LBFGSB
my_moop = MOOP(TR_LBFGSB)

# Massive 50-variable black-box optimization problem
# Completely hopeless for methods that rely on global models
for i in range(1, 51):
    my_moop.addDesign({'name': f"x{i}",
                       'des_type': "continuous",
                       'lb': 0.0, 'ub': 1.0})

# A simple convex simulation output with 2 outputs
def sim_func(x):
    xx = np.zeros(50)
    for i in range(50):
        xx[i] = x[f"x{i+1}"]
    # 25 variables that don't affect tradeoff, but need to be minimized
    tail = np.linalg.norm(xx[25:] - 0.5) ** 2 / 25
    # 25 variables that do affect tradeoff
    s1 = np.linalg.norm(xx[:25] - 0.2) ** 2 / 25 + tail
    s2 = np.linalg.norm(xx[:25] - 0.8) ** 2 / 25 + tail
    return np.array([s1, s2])

# Using a local surrogate to dodge the curse of dimensionality
# Notice that search_budget has to be greater than the number of variables
my_moop.addSimulation({'name': "MySim",
                       'm': 2,
                       'sim_func': sim_func,
                       'search': LatinHypercube,
                       'surrogate': LocalGaussRBF,
                       'hyperparams': {'search_budget': 200}})

# 2 objectives (using the single_sim_out library objective to minimize a
# single output of the simulation function)
my_moop.addObjective({'name': "f1",
                      'obj_func': single_sim_out(my_moop.getDesignType(),
                                                 my_moop.getSimulationType(),
                                                 ("MySim", 0))})
my_moop.addObjective({'name': "f2",
                      'obj_func': single_sim_out(my_moop.getDesignType(),
                                                 my_moop.getSimulationType(),
                                                 ("MySim", 1))})

# When solving big problems, it's often better to fix some acquisitions so
# we can focus on a few high-quality solutions
my_moop.addAcquisition({'acquisition': FixedWeights,
                        'hyperparams': {'weights': np.eye(2)[0]}})
my_moop.addAcquisition({'acquisition': FixedWeights,
                        'hyperparams': {'weights': np.eye(2)[1]}})
my_moop.addAcquisition({'acquisition': FixedWeights,
                        'hyperparams': {'weights': np.ones(2) / 2}})
# Keep one randomized acquisition to get some coverage of the Pareto front
my_moop.addAcquisition({'acquisition': RandomConstraint,
                        'hyperparams': {}})

# Turn on logging with timestamps
logging.basicConfig(level=logging.INFO,
                    format='%(asctime)s %(levelname)-8s %(message)s',
                    datefmt='%Y-%m-%d %H:%M:%S')

# 200 iterations * 4 acquisition funcs + 200 point search = 1000 eval budget
# This could take a few mins to run...
my_moop.solve(200)

# Display the values of x26, ..., x50 for all solution points
results = my_moop.getPF(format="pandas")
results[[f"x{i}" for i in range(26, 51)]].to_csv("local_method.csv")

# Plot results -- must have extra viz dependencies installed
from parmoo.viz import scatter
# The optional arg `output` exports directly to jpg instead of interactive mode
scatter(my_moop, output="jpeg")

The above code is able to approximately solve the problem on an extremely limited budget given the large dimension, and it produces the following figure of the nondominated points:

Scatter plot of solution points found for a high-dimensional problem

The solution is inexact, but the general shape of the Pareto front is already visible. Running for more iterations would further increase the accuracy.

Based on the problem definition in the code block above, a necessary but not sufficient condition for Pareto optimality would be x26=0.5, x27=0.5, … x50=0.5. To guess at our performance, the above method prints a csv file, selecting just columns x26, … x50 from the final dataframe:

,x26,x27,x28,x29,x30,x31,x32,x33,x34,x35,x36,x37,x38,x39,x40,x41,x42,x43,x44,x45,x46,x47,x48,x49,x50
0,0.5366946835438436,0.480535634334331,0.5160475188419098,0.6070197065831299,0.6672172711094266,0.48551619978110366,0.42610540131960145,0.4811010738722131,0.3425328910576567,0.7054110833892607,0.4352889180964643,0.7978816841760015,0.3063578175995894,0.5292602000450671,0.700800224549289,0.5377027107208002,0.25115698686137533,0.5288995346995834,0.4950069227982256,0.17614554775339386,0.43159599264801474,0.5490206585918878,0.6063607938295565,0.5025712376141976,0.47940143509202565
1,0.21201656695993382,0.3702045521967413,0.6190734162294054,0.6239469986512967,0.2862085782607271,0.517951134771002,0.5298746636054862,0.3585855775673113,0.3996404437543248,0.3506929113945627,0.5651786195889051,0.530138859271535,0.3864477482979602,0.6248489519955427,0.3162146560016008,0.24391249573503515,0.28071059514585706,0.370883066418232,0.5541698806411045,0.2265722859000449,0.47532448644242864,0.5388644174110008,0.6616590569894238,0.44813936868080606,0.4009390125097746
2,0.5561424758032281,0.44346734011216216,0.38998215633756317,0.4591350346310652,0.7328207110981155,0.44930979877679217,0.32144499533332593,0.39650142970701746,0.48863474948227376,0.507833114500235,0.26086792936556025,0.8056729955634335,0.33940034592710544,0.345840903491873,0.4053338723248634,0.3930238052957297,0.6466367621586258,0.5565063641229069,0.5887740374359774,0.4198769007664722,0.4827112230905737,0.24421871197214853,0.33210789218288583,0.39605677887110347,0.6878048361733072
3,0.44286300869997175,0.42778547451795623,0.35186496326814265,0.6513758449970273,0.4052515181173766,0.4078043109512294,0.5073185829990507,0.15106299313818577,0.4856048764844075,0.6413239003132395,0.6304134167473617,0.45121397858198714,0.5889145607293959,0.42748007301056445,0.40715934094130535,0.3592057831636401,0.4960409128034893,0.5380332043625424,0.36382354945134554,0.2013616120308228,0.49030851239578516,0.554521649474869,0.6327178428704407,0.41751222510187963,0.6772355123380417
4,0.23595077432304035,0.2732761441136762,0.6638340784257559,0.5120017574606377,0.2758737217174738,0.7126421327883947,0.275295911191266,0.7638257646563237,0.6382591391425461,0.48435491499904243,0.3626957409583479,0.5621126041659689,0.7893338288441889,0.52226623612655,0.3059973866695175,0.552189752202072,0.8135134811086226,0.567471337713878,0.5579479179061301,0.3894876873448983,0.3215337303042187,0.2700470849747288,0.5160923190693019,0.4965910930525222,0.5061571642517574
5,0.46149336815190034,0.629538330175257,0.7127560377416142,0.4442887226265191,0.4541224130287897,0.2633243966165434,0.24692939062063107,0.3784509135054398,0.26898224709879226,0.7357933931843939,0.3475382778291568,0.5171148295386534,0.3501570231861545,0.7028838956548591,0.716751049973538,0.5542627521702089,0.3166422001483687,0.732718408693968,0.44179341859052684,0.6915282387249422,0.45862007404657235,0.6557503406166848,0.7416320493350386,0.18700464597241073,0.4738179362452168
6,0.6492155188914149,0.7386180921623362,0.8191657705640728,0.6177120946296875,0.25545730090450314,0.6279436681962359,0.7456089531715299,0.4301192581459022,0.431919833147515,0.1825871560697172,0.6449108450032169,0.31818089007798733,0.595526081062856,0.5665665151006083,0.7741614283541227,0.36394801685054523,0.37142174861180277,0.270050607879592,0.30698220823439876,0.4400586749726125,0.7280106145094976,0.25166356859756034,0.7608913835265945,0.7186735352378809,0.21475152967771133
7,0.46511624903347193,0.5325359660189231,0.7363345308549895,0.47213948139686746,0.25221355324330574,0.702611596430013,0.19801360881641222,0.6751423733962472,0.5562234092469281,0.42163280717603746,0.7037384241492336,0.5266846511578003,0.35220590955802933,0.7865953485969641,0.5754496822029926,0.27427612361645637,0.6029163148214325,0.5187626047343203,0.8090796351007774,0.6298538773862838,0.3147845801056105,0.7979756538309042,0.33750508850227334,0.5539840566523027,0.28944326413397214
8,0.6240532503636728,0.8715552000151949,0.6698354610710465,0.3540450624148303,0.5379965591514205,0.523314143146974,0.17951857907235985,0.5644640557818029,0.7453920798128206,0.14683389329059543,0.6744902863826904,0.42902093888927334,0.6946287384636932,0.4387464460786994,0.6567829228190275,0.7714918137748795,0.20838260380861762,0.3866154762935723,0.7258152024877955,0.37958037618930507,0.4166675771106789,0.799775088460385,0.6424508553766978,0.5073058489695419,0.5122928196153836
9,0.8083228274693901,0.39564148332653515,0.35028907654874825,0.9668048714410297,0.2907483221250661,0.11999807062618467,0.5379576515967259,0.6673637446575738,0.2775098975155653,0.38931769193400007,0.534830206698696,0.37199641724695753,0.5451064223164794,0.1814521626676714,0.29179908599486015,0.7381192982516395,0.6697243534585621,0.5317930404167982,0.7379681323797941,0.64677195699539,0.7248551457955743,0.36593461099865,0.7152835594401916,0.41988317942116643,0.5214843561645818
10,0.4971316388237236,0.22682229991421982,0.5034779108983308,0.43868655362420184,0.7385335954546681,0.32204439723191597,0.9344507452250923,0.5779280376271245,0.6463574448323753,0.219346142470507,0.6288287743306861,0.3788840056203297,0.4041686697010589,0.4951429228829251,0.24735083357256366,0.4738042240031955,0.3974284988300453,0.5358793405580311,0.7142121969463976,0.8058424404118725,0.5819444460124938,0.541959646374663,0.45085970031475414,0.36754077273742913,0.3113884969912084
11,0.19951350287858372,0.7671634234797564,0.17989868355870525,0.34584115254781567,0.24481664558855418,0.3740867927207997,0.22618779201293282,0.40600617794191374,0.3074439894811378,0.5888098866502656,0.8233099259558825,0.6555812450709618,0.8236265332663156,0.3588332743741342,0.5425462065497779,0.26422902032702544,0.6440623209479265,0.5088102854045872,0.7009006177389732,0.5988893832521629,0.22836153221409056,0.48754020083308086,0.6534480841631513,0.2545829602727406,0.17811085551158462
12,0.619188400594239,0.6569994787428085,0.34386908022015206,0.3580530393730327,0.5519159178238913,0.7752589759093589,0.2882404268420529,0.5885724416310214,0.6484140329822388,0.25753817513149696,0.3786197761212532,0.21936649953049903,0.7588403580742802,0.573424933491924,0.49377388838693637,0.6049172970859933,0.29007401048194426,0.5191808323155452,0.5803150813931853,0.13360549678995692,0.7543773269703189,0.29244776310832327,0.16992855824317116,0.20649592811337603,0.7684992936901676
13,0.8089043454126528,0.24558629620503594,0.4156763954417,0.7750751509175771,0.5733692545340998,0.2679485326366635,0.53750309489663,0.33303474406978856,0.7084503201360435,0.8832895613058815,0.551111375044852,0.41344448359974495,0.401510817411759,0.3315761394811299,0.8119844473304381,0.49679670410916005,0.2331289666950867,0.3191770127712278,0.3905585632566304,0.8341357429713513,0.6363225606564503,0.7055674552769213,0.6318181997230172,0.6148893630916908,0.2473977765575856
14,0.3630333976468133,0.7312284036556909,0.4063879600787254,0.34682011182462447,0.43076419574801894,0.49614735963150935,0.5642599399841538,0.23592585554747167,0.4735232573773983,0.41591286744870337,0.7107062739910608,0.582879903830216,0.5784965960257064,0.6804472686056164,0.5271577265185332,0.42950454556205403,0.7094201476841995,0.15356846899708243,0.41489067248084943,0.6767144767496903,0.7291544282078588,0.3937279161140911,0.43346368619686954,0.6167654985691386,0.23988257785104986
15,0.6017643073886613,0.3922882389716225,0.6865410366409341,0.8524120136792644,0.33314740251996156,0.5517539198107381,0.2327922616398413,0.5536736125248068,0.5370815585082317,0.23074641605982843,0.435554111119259,0.5808195289088446,0.484419418046966,0.394608239223725,0.7211474951159841,0.375611436745661,0.7875268518268734,0.25454616319718004,0.45458907552995853,0.47339219099230107,0.4243546655769569,0.3260372408529185,0.6515562334746124,0.4616192647290729,0.5240954207844519
16,0.305441663990466,0.4256369598641956,0.5530375899730964,0.2833604667426151,0.51441960076015,0.5490441252909793,0.31899614978547386,0.39610768732306634,0.8387623919193579,0.5866002922680319,0.41109892150846616,0.7926501023816348,0.46490465122520463,0.5890330085147424,0.3337441558064368,0.8138882756594232,0.27921683389559543,0.4580925338430854,0.46433812752959863,0.4798901616824015,0.3568627007507853,0.27221780616665187,0.4295268033041352,0.6602843904712186,0.34315577280675214
17,0.48493200481135573,0.43609521301064386,0.3765575487986266,0.842615143516219,0.5137473362949972,0.7105427987004205,0.5959799744216097,0.5157955128852942,0.3549782073134401,0.40418583378308176,0.5831240472544407,0.4968322893613715,0.29510391175219086,0.7470204825699024,0.7292694773175746,0.49290601696963143,0.5269182911645355,0.5181521228053719,0.7265849547557288,0.6124109945880645,0.5432599236594075,0.33655171706368814,0.786333946818256,0.5276579600099066,0.6026212160010831
18,0.26073856777261273,0.34602946393585,0.5885191774152281,0.6017552570462076,0.36859034061765555,0.6226629674181852,0.36783675425109297,0.3522830393328373,0.4144792031204204,0.6621695010830078,0.7885543219711353,0.45723993940160584,0.31166586047143485,0.4538507250015391,0.3047875639565541,0.6051068870484023,0.24855959667438468,0.347584731181543,0.483278403568783,0.2745221614931336,0.717695929306571,0.3882041691579231,0.35282408842752033,0.7576750787772889,0.27371532305000945
19,0.6087577032918362,0.5163737338642169,0.7927612250769009,0.40107675259140624,0.25074158880469255,0.48682983665780055,0.7777176294352055,0.43527911526968577,0.7306670027820905,0.18286849588030185,0.4208555425616807,0.3320613413189294,0.5301087333811775,0.3562850384182146,0.47420942900080615,0.7022023824030751,0.23628440604611273,0.37753529090993965,0.615703375152078,0.5013486360338666,0.4763234219308978,0.6107173315558655,0.2209548846031878,0.5721368098253389,0.6748518330113652

Clearly, ParMOO has not found the exact solutions, but many of the solutions are quite close in all but a few columns of x26, …, x50. Whether these solutions would be accurate enough for a given application is entirely application dependent.

If you use these techniques in your research, consider citing our design paper, where we describe a similar example available in the solver_farm in Section 5:

@techreport{parmoo-design,
    author={Chang, Tyler H. and Wild, Stefan M.},
    title={Designing a Framework for Solving Multiobjective Simulation Optimization Problems},
    year={2023},
    institution={arXiv preprint},
    doi={10.48550/arXiv.2304.06881}
}