Re: DifferentialEvolution: custom mutation and recombination functions?
Karen Róbertsdóttir <[email protected]> Tue, 22 Aug 2023 04:31:30 +0000
| Newsgroups | gmane.comp.python.scientific.devel |
|---|---|
| Message-ID | <CAN=yiRyZ8T5ZBqnu7zmiyBR3L6hBa_sViwcisEq77A5tBzUFUw@mail.gmail.com> |
Any followup on this? Thanks! :) - kv, Karen þri., 15. ágú. 2023 kl. 21:30 skrifaði Karen Róbertsdóttir < [email protected]>: > Thanks for the response! Sorry for any misunderstandings I might have had > about the current internal architecture :) > > > The callable function would be completely > > responsible for generating trial vectors, i.e. doing the mutation AND > > recombination. > > No issues with that! > > > A possible call signature would be `strategy_func(candidate, population, > > rng=None)` > > Perfectly fine! But, question: > > > The strategy_func would be > > responsible for mutating (blending) members of the population together, > > doing the crossover/recombination itself, and returning a trial vector > with > > shape (N,). > > So it wouldn't be told what other candidate to perform recombination with > - it should pick recombination targets itself? I mean, that's workable, > just being clear on this. > > > Note that all entries in population are currently numbers in the range > [0, > > 1], and are scaled to 'actual values' using the lower and upper bounds. > If > > such a feature is added it would be reasonable for the strategy_func to > > receive population values in their actual range, i.e. scaled from [0, 1] > to > > [lowerlim, upperlim], sent to strategy_func, trial returned from func is > > then scaled back to [0, 1]. The to/from scaling would add some overhead. > > I can't speak for others, but I'm fine with receiving it in the [0, 1] > range and then scaling it myself, to avoid the need for the stock functions > to take that slight overhead hit. But whatever your preference is works > for me. > > > It's unclear what you meant by test functions returning a string of 1s .. > > There are no strings anywhere. > > Bad phrasing. "returning a np.array of 1.0s". I wasn't referring to a > literal string datatype. > > An example would be something like: > > ---------- > def strategy_func(candidate, population, rng=None): > return np.ones(len(candidate)) > > def minimization_function(x): > if np.all(x == 1.0): > print("Test passed") > sys.exit() > else: > return random.random() > ---------- > > ... which would be run with a low population size, ideally just one > candidate. So if it's working right, strategy_func returns a candidate > that's all 1.0, the minimization function sees it, and passes the test. > > > If such a feature was added it would be limited to generation of the > trial > > vector only. > > Correct. My apologies if my phrasing sounded otherwise. > > > It's not quite clear to me what your strategy function looks like, it'd > be > > interesting to see an example. > > Let me first describe an example from nature, before describing my case. > > The human eye has three different types of cones used for colour vision. > These are based on different opsin proteins, each of which has specific > genes responsible for their production. In particular, the M and L (green > and red) cone opsin genes are extremely similar, 96% identical. This > probably arose due to a tandem duplication - two copies of the gene arise, > one after the other. When this occurs, the initial impact is simply > increasing the amount of the original opsin protein production. But now > there's two separate genes, and they can drift independently without one > affecting the other. So now one drifts toward shorter wavelengths and the > other drifts toward longer wavelengths. Now the organism has the ability > to distinguish red from green, and can now see the difference between, say, > ripe and unripe fruit, or between a snake and a vine - conferring a > significant survival advantage! > > This couldn't have been possible without the tandem duplication event, > without one gene being able to turn into two. This cannot happen in scipy > as it exists now. But with a custom strategy_func, it could. If the coder > were representing "genes" in whatever their project is as clusters of > floating point numbers, they know how they're representing that data, and > can include gene duplication events, gene migration events, and so forth. > > Now, as for my specific example at present (note: this is not the first > time I've wanted such a feature - I use scipy.optimize.minimize a lot, I > love it! - but it's the first time I've bothered to get on the mailing > list :) ): > > One type of project I've used scipy.optimize.differential_evolution for > many times is developing 3d models for CFD (Computational Fluid Dynamics) > simulations - that is, to evolve optimal shapes for given tasks. Let's say > you wanted to evolve a wing - a candidate array might be of the format: > > [wingspan, chord, thickness, taper_ratio, twist_angle, incidence_angle, > dihedral_angle, sweep_angle .... ] > > You'd write a function to generate a 3d wing model from those set of > parameters, and then in your minimization function, you generate the model, > pass it off to e.g. OpenFOAM to run the simulation, save the forces from > the output, and then return a value based on those functions that > represents the performance of your wing. > > All well and good. And this is what I've done every time so far. But > there's a couple problems. The first is obviously that the mesh can only be > altered along the specific ways you design it to be altered - it can't > invent something innovative. And secondly, writing a function to generate > a mesh from parameters can be surprisingly time-consuming and challenging. > Including the fact that if the mesh happens to accidentally self-intersect, > it'll generate an aphysical model, and then the simulation can do all sorts > of crazy things. Every time you start a new project, you have to write a > new parameterized model-generation function from scratch. > > So this time around I decided to try something new: I want to make a > *generic* optimizer. Where there's no hard-coded model-generation function > at all for each task - where you can just provide an initial guess model, > and it can change it at will. Where one only has to provide constraints > and objectives for the simulation on which the results will be evaluated. > So then the question comes: how do you represent such a 3d model as a > candidate for evolution? > > A naive approach would be, "well, a mesh is vertices and faces, so let's > just list the vertices, then list their face indices, and call that good." > But more than a couple minutes thought shows that this is a disastrous > idea. Indices are integers. Mutating from one index to another is a > nonsensical change, and usually one that will make a broken mesh. And > crossbreeding/recombining meshes also makes no sense - even crossing > vertices will often yield incoherent results, let alone crossbreeding face > indices! > > Instead, I settled on an incremental generative approach. The first 14 > floating point parameters of the candidate describe an initial anchor face > for the model, and then each subsequent group of 11 floating point > parameters describes a new face to be incrementally added to an edge of the > model (each of these groups of 11 parameters can be thought of as a gene). > That is to say, each edge in the model has an edge anchor ID associated > with it, and all still-available edges on which the model can grow on are > in an unused_edge_id list. So for each gene (e.g. each new face), it > searches through unused_edge_id list, finds the one that most closely > matches its anchor ID, and builds itself there. To build a new face, it > grows a vector a fixed distance in-plane out of the midpoint of the edge > it's being attached to, and that vector is then rotated in-plane (e.g. > along the root face's normal) and out-of-plane (e.g. along the attachment > edge), to find where to add a new vertex. If the new vertex is close > enough to an existing vertex, it snaps to that vertex (allowing the model > to create closed shapes) - otherwise it adds a new vertex in place. A face > is then formed from the preexisting edge to the new vertex (or preexisting > vertex, in the case of snapping), any formerly unused edges are marked as > used, and any newly added edges have their anchor IDs added to the > unused_edge_id list. > > One of the 11 parameters for each gene is replication_count, so that a > model can evolve to make many copies of a given added face (each picking > the most similar remaining anchor ID, as per before). > > Once all the "genes" in the candidate have been processed, the mesh is > then run through bpy to use Blender to extrude it by a given thickness > (based on each face's material properties, which are also among the 11 > floating point parameters), and self-collision checks are also performed. > Any defective models return from the minimization function immediately > (returning an extremely high value, as a failure). Otherwise the model > then goes on to OpenFOAM for a CFD simulation as usual. > > (Current status: model generation seems nearly debugged, at least for > simple models. Haven't gotten into bpy extrusion or self-collision tests > yet. But I figured I should open this conversation re: scipy changes now > because they might take some time). > > Given the above, you can see how gene duplication (e.g. copying a cluster > of 11 floating point parameters and overwriting a different one) followed > by genetic drift has the potential for doing the same sort of thing that > happens in biological systems: allowing a piece of functionality to develop > into new pieces of functionality while simultaneously not destroying the > original functionality. If the replication_count on the gene for a given > facet description is 5, perhaps it overwrites some unimportant other gene > (genes expressing really tiny or unrealistically elongated facets late in > the genome could be preferentially targeted for overwriting) and leaves us > with one having a replication_count of 3 and the other with a > replication_count of 2 - aka, the same net result as before. But now these > two distinct genes can drift apart from each other and develop into new > types of functionality. > > (As a side note, Scipy's inability to save and resume the population > during differential_evolution optimization used to be really annoying, > given how long CFD optimization tasks take. However, I did find a kind of > cheap hack that I've been using ever since - since the random number > generator is deterministic, I simply have the minimization function create > a hash value for the candidate, and store the results of the simulation in > a hash table, which I save to disk. Then when I need to resume, I just > load up the hash table, and if a candidate has been encountered before, it > just immediately returns the previous run's simulation results rather than > re-running the simulation. It's an awkward hack, and wouldn't work on tasks > where the minimization function is really fast, but for slow tasks like > CFD, it works :) ) > > Custom strategy_funcs can be of course used for things that have nothing > to do with genes. For example: sometimes - as the docs note - a user may > want part of their candidates' data to be interpreted as integer data. How > do you mutate or crossbreed integers and have them make sense? Well, > that's really going to be task-dependent. Maybe the integer means "number > of iterations" - if so, then perhaps simple interpolation is best. But > maybe it's a category - in that case, interpolation is incoherent, and you > should either keep it the same or randomly pick a new category. And if it > is a category, and that category influences some other values in the > candidate, then that may affect how you want to alter those values. Maybe > if round(candidate[0]) is Category == 3, then you want the floating point > value at candidate[1] to be between 1.0 and 10.0, but if it's Category == 5 > then maybe you want candidate[1] to be between 1.0 and 5.0. Again, it's > task dependent. > > Honestly, there's no limits to what one could do with access to custom > strategy_func implementations. If one wanted, they could outright train a > neural network on candidate values and how well those candidates perform, > and let the neural network mutate candidates, so that the changes aren't > random, but are rather guided by complex statistics about "what sort of > alterations to a previously successful candidate are most likely to make an > even more successful candidate?". Think, say, a protein-folding > optimization task. > > (I have no plans to do such a thing personally, but it's just an example > of how far one could take this if needed for complicated, slow tasks) > > > > If there are no plans to implement it, I might (or might not, > depending) > > > be able to find the time to do so. Though I know 98% of the time > required > > > would not be coding / testing, but rather figuring out how to setup and > > > develop for the scipy test environment and figuring out how to > contribute > > > the changes > > > > THe time split would probably be 5% setup, 30 % writing code, 30% writing > > test cases, 25 % polishing. > > I like your optimism, and I'm sure it would be for someone who has > experience like you, but the last time I contributed to a project on > Github, just figuring out how to create the pull request for the finished > code took 40 minutes. ;) I've learned pessimism over the years over how > long it can take to set up dev environments and learn my way around them. > But maybe my pessimism is unjustified here. :) > > - kv, Karen > > _______________________________________________ SciPy-Dev mailing list -- [email protected] To unsubscribe send an email to [email protected] https://mail.python.org/mailman3/lists/scipy-dev.python.org/ Member address: [email protected]