New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Multiple solutions [feature proposal] #23

Closed
jeffreyhanson opened this Issue Sep 29, 2017 · 11 comments

Comments

Projects
None yet
2 participants
@jeffreyhanson
Copy link
Contributor

jeffreyhanson commented Sep 29, 2017

I've talked to a few people about the ability to output multiple solutions and it seems like a lot of people want it. Here are some ideas on methods to achieve this:

  1. Implement Bender's cuts (as in the raptr R package). Essentially, the problem is solved, and then additional constraints are added to forbid the solution, and the problem is solved again, and more constraints are added, and etc, until a desired number of solutions are obtained. The problem with this approach is that it is very computationally intensive with no shortcuts. To my knowledge, if the users requires n solutions the problem must be solved n times. In summary, good points: it can be used to obtain n solutions that are nearest to optimality and lets us use ROI (#16); bad points: very slow.

  2. Access the solution pool. Gurobi will store feasible solutions it comes across while solving a problem. These can be retrieved after solving the problem for little computational overhead. In summary, good points: quickly get multiple solutions; bad points: quality of solutions is unknown, requires writing C++ code to interface with Gurobi, and prevents us from using ROI to easily extend range of solvers (#16).

  3. Shuffle solutions prior to solving problems. It appears that simply re-ordering columns in an LP can cause different solutions to be returned. This is relatively easy to implement by randomly re-ordering columns prior to solving. In summary, good points: potentially obtain different solutions, all solutions meet the solution gap, plays nicely with ROI (#16); bad points: no guarantee solutions will be different, very slow as all problems must be pre-solved and then solved.

3a. Shuffle solutions prior to solving problems but after pre-solving problems. Gurobi seems to spend most of the time during the pre-solve step, so avoiding having to repeat this step is most desirable when generating multiple solutions. One potential solution might be to use custom C++ code to (1) presolve the problem, (2) shuffle the columns, (3) solve the shuffled problem, and repeat steps 2-3 n times for n solutions. In summary, good points: all points in 3 but much faster than 3; bad points: idea remains to be tested, will not play nicely with ROI (#16).

Thoughts? Does anyone have any other ideas, or preferences for particular approaches? I suppose with regards to ROI in 2 and 3a, we could just implement a special case for Gurobi but I think the main benefit of using ROI is that we don't have to interface with solvers directly.

@ricschuster

This comment has been minimized.

Copy link
Member

ricschuster commented Sep 29, 2017

Hi Jeff,

This is an important point and outputting multiple solutions is certainly desirable. I personally prefer 1, even though its computationally intensive, because the way we have generated multiple solutions for far was to vary targets or constraints to create a pool of options and solutions, before looking at the overlap between them.

The rational behind this is that we assume that the ILP solver finds the optimal solution to a problem, so there is little reason to create a solution pool from the same problem, as the solutions in the pool will most likely be very similar. This might be a misunderstanding on my part, but that's my assumption.

The way that I have approached 1, at least via Shiny in the past was to allow users to either enter parameters on targets and/or constraints manually or upload a csv file where the scenarios are specified. Example here: http://forbasin.forestry.ubc.ca/CDFCP.v0.22.3/
(I think currently the upload doesn't work, so you will have to create examples manually)

This is also the way I was planning to implement this with the prshiny function. I think for practitioners this way would make the most sense, but I haven't explored 2 or 3, so am not sure if that could be feasible as well.

@jeffreyhanson

This comment has been minimized.

Copy link
Contributor

jeffreyhanson commented Sep 30, 2017

Hmm, you raise some really good points. I think the best strategy probably depends on why the user wants multiple solutions.

If the user wants the top n solutions nearest to optimality and views optimality is being really important, then option 1 is most useful for them. However, for medium sized problems the top 100 solutions will probably be very similar with only a few planning units that are different between them, so I don't know how useful this will be in practice?

My impression is that many users would like a portfolio of solutions that are very different in terms of the planning units selected, but still within some specified level of optimality. It seems to me that a lot of conservation planning isn't actually done with "real" cost data---though I don't have any experience working in real-world conservation planning --so optimality isn't really considered too important so long as the solutions aren't really poor.

What do you think? I could try out 3 or 3a with some simulated data and see how it goes.

@ricschuster

This comment has been minimized.

Copy link
Member

ricschuster commented Sep 30, 2017

This is a tricky one. Maybe we should ask the Marxan mailing list for input on this?

To me the main difference between 1 and 2/3 is, that in the case of one the scenario generation is up to the user, in 2/3 its more similar to Marxan in the sense that the user specifies one problem that is tried to be solved a number of ways by the solver.

The main issue I have with 2/3 is that we would have a way to get to the optimal solution with ILP, but we deliberately choose to forgo reaching an optimal solution for the sake of generating a portfolio of solutions. To me that seems a bit strange, but I'm not sure how others think about this. I also wonder how different solutions would be if we say choose a pool of solutions that have a 5% gap. That gap might be acceptable for most (could also let the user specify the gap), but we would need to test this to see how much variation we actually get in the pu selections.

It also sounds like that for both 2 and 3 we would have issues using ROI, which is a shame, because it sounds like that ROI would be good to implement moving forward.

Maybe try 3 or 3a, see how it goes, and also ask the mailing list?
1 is computationally intensive, but would require little effort in terms of code.

@jeffreyhanson

This comment has been minimized.

Copy link
Contributor

jeffreyhanson commented Oct 2, 2017

Yeah, asking the Marxan emailing list sounds like a great idea.

I think 1 & 3 are equally computationally intensive. They both involve resolving the entire problem n times. Both 1 and 3 are compatible with ROI. It's just 3a that could potentially offer computational shortcuts and break compatibility with ROI.

@jeffreyhanson

This comment has been minimized.

Copy link
Contributor

jeffreyhanson commented Oct 2, 2017

I've implemented option 3 in the multiple_solutions_3 branch. I've been exploring it's functionality, and have put together a minimal rmarkdown doc (ms.zip). I'd be really keen to hear what you think.

I also realised that option 3 could potentially be much faster than option 1 because each solution can be solved in parallel, which is useful because Gurobi spends most of its time in the presolve stage for simple problems (eg. classic marxan problems) does not use multiple threads (even if specified in to use multiple threads).

If you think this is a viable option for generating multiple solutions, I would be interested to hear your thoughts on implementing it. For example, should users generate solutions by running solve multiple times? Or should there be a number_solutions parameter in add_gurobi_solver? How would you see this interfacing with #22 and #25?

@ricschuster

This comment has been minimized.

Copy link
Member

ricschuster commented Oct 2, 2017

This is awesome Jeff! I've played around the the markdown doc you put together a bit and I like the way you set this up. I think its very easy to understand for prioritizr users to use this. I have added a few things to the markdown doc myself (MS_3.zip).

The first thing I did was to quickly put together a figure showing the pu selection frequency of the 1000 runs you had in the example. There are some options in the solution space for sure, which is I think a good starting point. I did also go ahead and reduced the gap to 5%, just to see how that would play out. Pretty decent I think.

One thing to note is that there is a lot of repeats in the solution space at the moment (at least for 5%). The question is, if those are true repeats or if the selection pattern just happens to be the same. Do you think there is a difference between those two or does it really mean that each time the exact same pu's are selected the solution is the same? I'm just wondering if maybe getting to one solution could be done via different routes and if that would still mean the solutions are the same or just happen to select the same pu's? If they are really the same, do you think there is a way to reduce computational load by removing duplicate setups before they are actually solved?

The speedup by solving in parallel sounds very promising to me.

I haven't tried this with one of my bigger datasets yet, but if you think the parallel processing could be implemented quickly, it would be great to give that a try with a bigger dataset.

I think a number_solutions parameter in add_solver would make the most sense. This way the user can use one function to specify the gap and the number of solutions they want to get out.

#22: I wonder if we should recreate the output structure from Marxan, where we export each solution and the selection frequency? This might get a bit crazy for bigger problems though. An alternative would be to have the selection frequency as a spatial output if inputs are spatial and the solution matrix as a csv. To link to the raster we could also create an index raster that gets id values which are used as the first column in the solution matrix csv file. In some ways I'm thinking the Marxan output is a bit outdated, but if we mirror Marxan's output structure here, users would have an easier time switching?

#25: A list of attributes for the solution matrix? I guess if we already have a solution matrix as output (I don't think a solution raster stack makes the most sense here), we could also create an output class with slots for 1) spatial representation, 2) solution matrix, 3...n) attributes from #25. What do you think?

I do like this solution space generating a lot, we just have to make sure this is sensible and makes computational sense. Its also great in the sense that one can get the optimal solution fairly quickly and can then relax the gap. This is quite different from Marxan and a huge improvement I think.

@jeffreyhanson

This comment has been minimized.

Copy link
Contributor

jeffreyhanson commented Oct 2, 2017

Yeah, the random shuffle method (3) does return a lot of duplicate solutions. I guess there's no real way to know how well the (3) is sampling the range of solutions within the gap without actually implementing Benders cuts (1) to determine the number of solutions within the gap. In other words, 300 unique solutions might sound good, but it's hard to know if this is 40 % or 100 % of the solutions within the 5 % gap.

Yeah, it's entirely possible that the duplicates could be returned via different algorithms/methods. But I don't think this matters too much. If Gurobi returns two identical solutions, they're still identical regards of whether it used algorithm A or algorithm B, so I would still remove all duplicates when calculating things like selection frequency. Yeah we might be able to re-use variable states to reduce incidences of duplicates, but then that means that each Gurobi run is dependent on previous runs, so this would make things tricky when running multiple Gurobi runs in parallel to generate the portfolio.

Yeah I think adding it to the add_*_solver makes the most sense. If that's the case, maybe we could add options (1) and (3) and use an additional argument to specify the method for generating multiple solution? So the function definition might look something like this?

add_gurobi_solver <- function(x, gap=0.1, time_limit=.Machine$integer.max,
                                                presolve=2, threads=1, first_feasible=0, number_solutions = 1L,
                                                multiple_solution_method = c("shuffle", "cuts"), 
                                                remove_duplicate_solutions = TRUE)) {....}

Re #22: Yeah I think the Marxan output probably getting a bit dated, but I think it's really important to make this identical to Marxan output. I think we risk annoying a lot of users if the outputs were different in any way.

Re #25: Yeah we could do that. Although I wonder if we could store additional solutions in extra columns in the data.frame or Spatial*DataFrame or as extra layers in a RasterStack. This would make it much easier for the users to access multiple solutions, than if we kept the extra solutions in an attribute. If we did that, then I think we would only need to store the properties for each solution (as a separate vectors?) in the attributes: objective values, solver return codes, run time, solution method?

@jeffreyhanson

This comment has been minimized.

Copy link
Contributor

jeffreyhanson commented Oct 2, 2017

I just had another thought. Instead of adding arguments to the solver functions, what if we added functions like add_cuts_portfolio(...) or add_shuffle_portfolio(...) that users could add to ConservationProblem class objects to specify how multiple solutions could be added?

This way users could generate multiple solutions using something like this:

  p <- problem(sim_pu_raster, sim_features) %>%
    add_min_set_objective() %>%
    add_absolute_targets(targ) %>%
    add_binary_decisions() %>%
    add_gurobi_solver(gap = 0.2) %>%
    add_cuts_portfolio(number_solutions = 10)
@jeffreyhanson

This comment has been minimized.

Copy link
Contributor

jeffreyhanson commented Oct 2, 2017

I think this has several benefits: it separates the portfolio code from the solver-specific code (which will be easier to maintain and make porting to ROI easier), it removes redundant combinations of arguments (eg. remove_duplicate_solutions and Bender's cuts), and also means that it will be easy to add new strategies for generating multiple solutions if we can think of any more. What do you think?

@ricschuster

This comment has been minimized.

Copy link
Member

ricschuster commented Oct 3, 2017

re #25: I agree, with your suggestion and might have been unclear before. It might be getting users into trouble for very large >1M pu's and >1000 features to have a raster stack comprised of 100+ rasters, but I would assume that users with that kind of size problem would know how to deal with problems like that.

I like the add_cuts_portfolio(...) suggestion as it does help to keep things cleaner and make adding new strategies easier.

@jeffreyhanson

This comment has been minimized.

Copy link
Contributor

jeffreyhanson commented Oct 3, 2017

Yeah, I agree.

Ok awesome. I've started implementing this along with #25, so hopefully this functionality will be in prioritizr soon.

jeffreyhanson added a commit that referenced this issue Oct 11, 2017

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment