Thanks for posting this. I just gave it a skim through, which makes me think I probably won't have much intelligent to say until the user docs are posted. I'll probably come back and read through again at that point.

That's the thing--it's such a big PR that it's hard to know what to comment on. I spent a while looking into ways I could maybe break this up into a few smaller changes, but it just wasn't practical, ultimately.

If any part of this should have been separated out, it's f456a8a; this PR very much depends on it, but that change could have been made first. It was originally part of another branch where there's more work being done to simplify how parameters are handled (in such a way that will also, finally, support Quantity parameters). But I haven't had time to get that 100% working. The part pulled out in f456a8a was very helpful to have here though, more than anything else.

@nden's fix in 0d66a05 addresses some of the known bugs, though there is another bug involving compound models using AffineTransformation2D that I need to address.

@embray I've noticed that now a sliced model does not have the correct parameters (perhaps a regression?)

offx = Shift(1, name='x')
offy = Shift(1, name='y')
sclx = Scale(1, name='xscale')
scly = Scale(2, name='yscale')
rot = Rotation2D(angle=60, name='rotation')
m = (offx & offy) | rot | (sclx & scly)
sm=m['xscale': 'yscale']
sm.parameters
Out[61]: array([ 1.        ,  1.        ,  1.04719755,  1.        ,  2.        ])

Evaluation of sm fails.

That's odd--I'm sure that used to work, so it might be a regression, though I'm not sure where it was introduced. Probably when I reworked how compound models of model instances are handled.

I haven't quite gotten to the bottom of this yet but I think it's a minor error somewhere. I have to head out for most of the rest of the day, but I think I'll be able to fix this tomorrow.

Another interesting thing that doesn't work with this example is:

sm = m[:'xscale']

in other words, a slice using every part of that model except for the final yscale. That example was never meant to work--the end result is something like:

(Shift & Shift) | Rotation2D | Scale

where the result of Rotation2D, which has two outputs, is being passed to Scale which has only one input, so it doesn't make sense to begin with and is expected not to work. There are two possibilities for how to handle this:

1. When taking such a nonsensical slice, simply raise an exception immediately.
2. It might be possible (and I think nicer) to implicitly convert this to the equivalent of:

(Shift & Shift) | Rotation2D | (Scale & Identity(1))

So the x output from Rotation2D is scaled, while the y output is passed through. I'm not sure how easy that would be but it's worth looking into. In the short term making it erroneous (option 1) will probably be easier.

Yeah, I was also thinking how to handle these cases. Another one is sm = ['y':'rotation'] . There are other ways for this to fail depending on the type of model combination, so may be requiring users to be explicit and raising an Exception is better than trying to guess the correct thing and supply an Identity model. It's a matter of matching n_inputs with n_outputs, right?

Yes--currently when combining two models in an expression there actually is no explicit checking on input/output compatibility which obviously should be. That's one of the items on my TODO list still. Since that's pretty easy I'll probably put it in before merging--should just take a few minutes. I agree that explicit is better than implicit, though at the same time I think putting in an Identity is fairly straightforward as well and might seem natural to some. I'm not sure though--if someone else thinks that's a good idea it can be changed later.

Is this for 1.0 or later?

It should be for 1.0.

@cdeil - this is one of the high priority features for 1.0 (see https://trello.com/b/XzwPX49N/brainstorming-for-1-0). It would be nice if it could be included, if the current issues can be solved.

@embray - I see you are working on the documentation now. I can review and test this tomorrow morning with several examples once the docs are up.

Weird...it seems like some of the doctests ran on Python 3? I thought we didn't enable that yet... I don't have those failures locally.

Tests in docstrings do seem to be running in Python 3 by default (saw this earlier this evening on a different PR, locally). Tests in the docs don't run yet.

I've tracked down the main bug in slicing. It's extremely subtle. Updating coming soon.

Meanwhile it seems like there's a bug with my isinstancemethod function on Python 3, which could be the root of some other issues. It's strange though, I tested that very extensively on all Python versions when I first wrote it so I don't know what's causing it to fail now.

There's still the [doctest] astropy.utils.introspection.isinstancemethod failure on Python 3:
https://travis-ci.org/astropy/astropy/jobs/44530894#L2584

And the sphinx build fails because of these warnings (see https://travis-ci.org/astropy/astropy/jobs/44530890):

/home/travis/build/astropy/astropy/docs/api/astropy.modeling.Model.rst:131: WARNING: error while formatting arguments for astropy.modeling.Model.rename: <unbound method _ModelMeta.rename> is not a Python function
/home/travis/build/astropy/astropy/docs/api/astropy.modeling.Model.rst:175: WARNING: error while formatting arguments for astropy.modeling.Model.rename: <unbound method _ModelMeta.rename> is not a Python function
looking for now-outdated files... none found
/home/travis/build/astropy/astropy/docs/api/astropy.modeling.Parameter.rst:12: WARNING: py:obj reference target not found: Parameter.__init__

This is huge ... thanks @embray! I will try this out today.

Looks like composite models always make copies of their parts.

Assume you have a simple signal model with a few parameter and a complex background model with dozens of parameters:

signal = SignalModel()
background = BackgroundModel()
compound = signal + background
# fitting
compound._tree.left.value # want to access signal model only.

Is there a way to do this without accessing the private _tree member and maybe by not having to navigate the tree myself, e.g. by model part name?
Also printing a string representation of the tree (in this case signal + background) would be very nice for complex composite models.

One thing I'd be very interested in is how to write a composite model in a library (e.g. in Astropy itself or an affiliated package) and have the parameter names and docstring turn out OK?

I'll try to do this tonight, but if you could put an example in the docs how to do this, that would be great!

@cdeil

There's still the [doctest] astropy.utils.introspection.isinstancemethod failure on Python 3

Yeah, I will fix that too. I was surprised that those doctests were running on Python 3 though--IIRC we hadn't enabled that yet, but I could be mistaken (maybe it was just the .rst doctests we didn't enable yet). Anyways, that one is failing because it uses __metaclass__ which won't work on Python 3, so I need to use six.add_metaclass for that example, I guess. Not ideal but oh well.

I'll fix the other doc issues as I work on updating the docs for this.

Looks like composite models always make copies of their parts.

That's correct. If needed, we can later look into ways to maybe get around this, but for now it's the safest default choice.

Is there a way to do this without accessing the private _tree member and maybe by not having to navigate the tree myself, e.g. by model part name?

That is a good question, and one I've brought up myself a few times (see for example here) There are a few possibilities right now. First of all, I'll note that when you make a compound model out of classes, the class __repr__ shows some useful information. First let me say something about classes, then I'll get back to your specific example. Using your example, it would look something like:

>>> Compound = SignalModel + BackgroundModel
>>> Compound
<class '__main__.CompoundModel53'>
Name: CompoundModel53
Inputs: ('x',)
Outputs: ('y',)
Fittable parameters: ('param_0_0', 'param_1_0', 'param_0_1', 'param_1_1')
Expression: [0] + [1]
Components: 
    [0]: <class 'SignalModel'>
    Name: SignalModel
    Inputs: ('x',)
    Outputs: ('y',)
    Fittable parameters: ('param_0', 'param_1')

    [1]: <class 'BackgroundModel'>
    Name: BackgroundModel
    Inputs: ('x',)
    Outputs: ('y',)
    Fittable parameters: ('param_0', 'param_1')

You can see a few things going on here--I think most of it is clear, but I am definitely open to better suggestions for how to format this output. Right now it's just what was useful enough for debugging purposes, but I think there are all kinds of possibilities for something slicker (especially in the IPython Notebook). You can also see that the compound model class has a pretty uninteresting generated name. We might be able to make more interesting default names, though that could also get a bit out of hand for complicated expressions. However, you can also give your compound model a custom name like so:

SignalPlusNoise = (SignalModel + BackgroundModel).rename('SignalPlusNoise')

or another, I think sort of nice syntax is:

class SignalPlusNoise(SignalModel + BackgroundModel):
    """Add a docstring too, why not?"""

As for navigating, once you've printed out the class and see the expression and its components you can access the individual components by their index:

>>> Compound[1]
<class 'BackgroundModel'>
Name: BackgroundModel
Inputs: ('x',)
Outputs: ('y',)
Fittable parameters: ('param_0', 'param_1')

or if the components have custom names, you can use the names instead of the numerical index, which is of course usually preferable and a good reason to give custom names (or sometimes I refer to them as labels).

All of the above applies when combining model instances like in your example:

>>> signal = SignalModel()
>>> background = BackgroundModel()
>>> compound = signal + background

The only problem with that right now is that the __repr__ doesn't display as much useful information. There's still a compound model class underlying compound, so you can get that info by entering:

>>> type(compound)

but I tend to agree that compound model instances should display more useful information in their own right. As for navigating the expression tree, as I wrote in the comment I linked to, I'd be open to ideas for how to better expose that (I think for most cases looking at the __repr__ and then using indexing is pretty good, but that still doesn't give a good API for introspecting the expression itself).

One thing I'd be very interested in is how to write a composite model in a library (e.g. in Astropy itself or an affiliated package) and have the parameter names and docstring turn out OK?

Right now there's not enough in the way of controlling parameter names. That's another item that's high on the TODO list, but was a bit too tricky to resolve for the first version. There are other ideas floating around associated with this. @ChrisBeaumont had the idea, which I like, of also allowing semantic hints on parameters. For example, many of the models have an "amplitude"-like scaling factor--if multiplying multiple models together it might make sense to combine their amplitudes into a single parameter, and such semantic hints would allow that. Other categories of parameters he identified included "width", "offset" and "slope", and probably others. I think it's an idea worth pursuing, though the details are still a bit hazy.

In the meantime I think providing custom labels for parameters on compound models is not far out of reach. I'm open to ideas for that.

Also, in the meantime, when you combine multiple models it generates some parameter names based on the names of the parameters in the sub-models, along with an index attached to it. This is a bit of a "good enough for now" kludge until there's something better. Though it's also possible to access parameters a bit more nicely by drilling down using indexing. For example if model A has a parameter amplitude and we write:

>>> m = A(1, name='A') + B() * C()

then we can access that parameter in two ways:

>>> m.amplitude_0
Parameter('amplitude_0', value=1.0)
>>> m['A'].amplitude
Parameter('amplitude', value=1.0)

These are just two separate "views", if you like, of the same parameter (updating one will update the other and vice versa). Ideas for how to make this even better are welcome.

If there are no further comments on this I'd like to merge soon-like so we can have it for AAS?

If there are any followup comments once more people get a chance to look at it we can address those later I think. Not all of this is set in stone.

@embray - yes, please feel free to merge this, as I think it is definitely ready for user testing. Thanks for all your work on this! 😄

🎉