Skip to content
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

Add support for refitting and nuisance averaging #343

Closed
wants to merge 7 commits into from

Conversation

kbattocchi
Copy link
Collaborator

No description provided.

Copy link
Contributor

@heimengqi heimengqi left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Mostly looks good to me. But I all agree on what you've mentioned in our call, it's annoying to loop all possible arguments and add a setter of it, and it's hard to maintain. We might need to brainstorming and see is there any alternatives.

econml/_rlearner.py Outdated Show resolved Hide resolved

@property
def discrete_treatment(self):
return self._discrete_treatment
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we really need to expose the first stage related arguments and allow users to set something else? I think the value of refit is to save computation time when the nuisance models are the same and the data and results are cached. Even they reset a value here they can't only refit final stage but fit the entire pipeline again, it's the same as reinitialize a new estimator and call fit in term of computing time and lines of code. Also it's wired to allow user to change the discrete_treatment flag since T is fixed and it should be either continuous our categorical. If we only allow user to change final stage specific argument, it could also simplify tons of properties and setters across multiple classes?

Same for _RLearner and classes under DML file.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

My idea was to make all things that you can set to some value in the initializer also available to update at any other time. It's true that many of these are not necessary for the refit functionality alone, and also that a workaround would be to create a new instance.

However, from the user perspective I think sklearn's behavior (where you can set any attribute at any time) would be much easier to understand than a policy where you can only set the attributes that don't affect the nuisances, since it is not always completely obvious which attributes affect the nuisances (e.g. the featurizer affects the Y residuals if we use the linear_first_stages flag, but doesn't affect any nuisances if it's False).

Eventually, I think we should enable setting attributes even for estimators that don't have support refit (e.g. the metalearners), but for this PR I'm just starting with the OrthoLearner class hierarchy.

def model_final(self, model):
raise AttributeError("LinearDML final model can't be chnaged from "
"StatsModelsLinearRegression(fit_intercept=False)")

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We can't change final stage model but we could change some arguments under that model? e.g. cov_type, like what you did for SubsampledHonestForest

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Another typo here I just randomly saw it, "changed"

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, I think that would make sense - but for now, cov_type is set based on the inference object's cov_type, it's never set directly by the user, so I just left that behavior as is for now. Ultimately this ties back into our discussion of whether to change the DML hierarchy by unifying the LinearDML and SparseDML classes and cleaning up how we implement inference, but that's out of scope for now.

return self._discrete_treatment

@discrete_treatment.setter
def discrete_treatment(self, discrete_treatment):
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Does this need to trigger a changed in the discrete treatment transformer??

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

As discussed offline, I don't think so because the transformation happens at fit time.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I found a related bug, not sure if related to the transformer. So the following works:

est = LinearDML(model_y=RandomForestRegressor(),
                model_t=RandomForestClassifier(min_samples_leaf=10),
                discrete_treatment=False,
                linear_first_stages=False,
                n_splits=6)
est.fit(Y, T, X=X, W=W)
te_pred = est.effect(X_test)
te_pred[0]

But the following raises an error, while from the current semantics it's supposed to work fine:

est = LinearDML(model_y=RandomForestRegressor(),
                model_t=RandomForestClassifier(min_samples_leaf=10),
                discrete_treatment=True,
                linear_first_stages=False,
                n_splits=6)
est.fit(Y, T, X=X, W=W)
est.discrete_treatment = False
est.fit(Y, T, X=X, W=W)
te_pred[0]

Which means that there is some state of the object that is maintained but should have been reset when resetting the discrete_treatment.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hm. I wouldn't expect that to work since it uses a classifier for model_t but a continuous treatment. I think changing discrete_treatment would require that either the model_t is 'auto' (although I need to fix this logic) or else that the model also changes.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thats why I added the first example that when called a priori it works. The classifier can handle the discrete treatment. So there seems a problem with changing diwcrete treatment

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Turns out this was due to a typo in property name for the DML setter. It will work after my next push.


@_OrthoLearner.discrete_treatment.setter
def discrete_treatement(self, discrete_treatement):
# super().discrete_treatment = discrete_treatment
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In general: we shoulld change the name of the _RLearner.model_final to _RLearner.rlearner_model_final, and _OrthoLearner.model_final to _OrthoLearner.ortholearner_model_final and then we should just be setting these variables. That would clean up the code even if the refit stuff were not there.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Agreed.

inference=inference)
cache_values=cache_values, monte_carlo_iterations=None, inference=inference)

@DML.model_final.setter
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I feel that for all these setters, there should be a logic within a setattr, as opposed to a separate method for each one, just so that we reduce lines of code. But this is secondary and I'm ok with doing it post release.

@@ -1065,3 +1064,59 @@ def test_deprecation(self):

d = pickle.dumps(LinearDMLCateEstimator())
e = pickle.loads(d)

def test_refit(self):
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think we need several more tests here that the new params affect the final stage. Things like changing fit_intercept and then testing that all our effect methods take the change after refitting.

Similarly maybe some other change in sparse_linear (e.g. alpha) or n_estimators in ForestDML.

Also a test that the "score_" attribute changes accordingly to the new model. Maybe do a model_selection experiment, where you run refit with two-three hyperparams, with 1 of them being the obvious best and the rest leading to random noise and read the score_ method each time after fit, and then make sure that the correct final model is chosen.

Similarly, that we can also do "score(Y, T, X, W)" out of sample and that the score returned is based on the new model.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Also maybe a test:
If I am in SparseLinearDML and I do:
SparseLinearDML().model_final = LinearRegression()

would that make the final model be LinearRegression()?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Finall, I think we definitely need a "check_is_fitted" to be called in all our predict methods and raise an error if not fitted. This is even more important now, as when I change a param we need to be invalidating the fit and so if someone calls effect, they shouldn't be getting the old fitted model final results, but an error that says that estimator is not fitted yet.

Currently our error messages are random internal errors of the estimators that some attribute is missing and not an interpretable error of the form "estimator is not fitted yet!".

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think we should use an approach similar to sklearn, where fit sets internal attributes of the fit that are stateful attributes and doesn't alter the inptu params. (e.g. at fit time we create a model_final_ which is the fitted model final, we also set models_y_ which are the fitted models_y etc, these are all attributes created by fit and the inptu params remain unaltered as static values).

Then sklearn has a generic check_is_fitted: https://scikit-learn.org/stable/modules/generated/sklearn.utils.validation.check_is_fitted.html
which just checks that the object has attributes ending in underscore.

For now we could go with simpler solutions by just setting a flag: "fitted" etc.

Copy link
Collaborator

@vsyrgkanis vsyrgkanis Dec 20, 2020

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Post-release we also definitely need to showcase this functionality in the dml notebook. For instance having a cell that does sth like the following as it is a superb functionality!

from econml.sklearn_extensions.linear_model import StatsModelsLinearRegression, StatsModelsRLM, DebiasedLasso
from sklearn.linear_model import ElasticNetCV, LassoCV
from sklearn.preprocessing import PolynomialFeatures

est = DML(model_y=RandomForestRegressor(min_samples_leaf=20, max_depth=3),
          model_t=RandomForestRegressor(min_samples_leaf=20, max_depth=3),
          model_final=StatsModelsLinearRegression(fit_intercept=False),
          random_state=123)
est.fit(Y_train, T_train, X=X_train, W=W_train, cache_values=True)

print(list(zip(est.coef_, *est.coef__interval())))
print(est.score_)

est.model_final = StatsModelsRLM(fit_intercept=True)
est.refit()
print(list(zip(est.coef_, *est.coef__interval())))
print(est.score_)

est.model_final = DebiasedLasso()
est.refit()
print(list(zip(est.coef_, *est.coef__interval())))
print(est.score_)

est.model_final = Lasso(alpha=0.1)
est.refit()
print(est.coef_)
print(est.score_)

est.model_final = ElasticNetCV(cv=3)
est.refit()
print(est.coef_)
print(est.score_)

est.featurizer = PolynomialFeatures(degree=3)
est.refit()
print(est.coef_)
print(est.score_)

I wish we could also use seamlessly the NonParamDML variants and the CausalForestDML variant, but those have a different rlearner_model_final, so it seems harder. but that would be very useful.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe if we just had an input flag that says: linear_model_final, then if that is True, we just don't use the weight trick and all functionality as in DML is good to go. and if linear_model_final=False, then we use the weight trick and we raise an error in the case when the user passes multiple treatments, that this is not allowed (I blieve that this will actually be done automatically but the model final wrapper).

So I think it's just a matter of adding the flag: linear_model_final, and then we could do:

est.linear_model_final=False
est.model_final = RandomForestRegressor()
est.refit()

This would allow us to do everything through refitting (except CausalForestDML, whcih uses a separate model final wrapper).

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For now we can do some of this by using the NonParamDML, but this still requires that the model accept weights. While for linear models we don't need weights. That seems good enough for now, but maybe the linear_model_final flag is still a good idea.

from econml.dml import NonParamDML
from econml.sklearn_extensions.linear_model import StatsModelsLinearRegression, StatsModelsRLM, DebiasedLasso
from econml.sklearn_extensions.ensemble import SubsampledHonestForest
from sklearn.linear_model import LassoCV
from sklearn.ensemble import RandomForestRegressor, GradientBoostingRegressor

est = NonParamDML(model_y=RandomForestRegressor(min_samples_leaf=20, max_depth=3),
                  model_t=RandomForestRegressor(min_samples_leaf=20, max_depth=3),
                  model_final=SubsampledHonestForest(),
                  random_state=123)
est.fit(Y_train, T_train, X=X_train, W=W_train, cache_values=True)

print(list(zip(est.effect(X_val[:1]), *est.effect_interval(X_val[:1]))))
print(est.score_)
est.model_final = DebiasedLasso()
est.refit()
print(list(zip(est.effect(X_val[:1]), *est.effect_interval(X_val[:1]))))
print(est.score_)

est.model_final = Lasso(alpha=0.01)
est.refit()
print(est.effect(X_val[:1]))
print(est.score_)

est.featurizer = PolynomialFeatures(degree=3)
est.refit()
print(est.effect(X_val[:1]))
print(est.score_)

est.featurizer = None
est.model_final = RandomForestRegressor(min_samples_leaf=20, max_depth=3)
est.refit()
print(est.effect(X_val[:1]))
print(est.score_)

est.model_final = GradientBoostingRegressor(min_samples_leaf=20, max_depth=3)
est.refit()
print(est.effect(X_val[:1]))
print(est.score_)

Copy link
Collaborator Author

@kbattocchi kbattocchi Dec 21, 2020

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Regarding checking fittedness, I agree that we should make sure that fit has been called before calling effect (and we don't always do this currently). However, note that sklearn will not reset an estimator to "unfitted" if attributes are set in between calls to fit; for example:

  1. fit an estimator
  2. set an attribute
  3. call predict

will work just fine, and generate predictions using the previously fit model (using non-current settings). I'd expect our behavior to match.

econml/_ortho_learner.py Outdated Show resolved Hide resolved
self._cache_invalid_message = None
return self

# TODO: this doesn't currently allow inference, because that is handled by the wrap_fit attribute
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So is inference the same object as the one called by fit?

I guess this is also related to the fact that the inference object does not have the new model final. I suspect this might need to be solved simultaneously (i.e. propagate model final to the inference object and allow refit to take as input "inference=...")

fitted_inds = None

for _ in range(monte_carlo_iterations or 1):
nuisances, new_inds = self._fit_nuisances(Y, T, X, W, Z, sample_weight=sample_weight, groups=groups)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If we fit multiple times, what version is models_nuisance storing? Is it storing for each mc iteration and for each fold in the kfold, the fitted model? Or is it storing only the final version of mc iterations? The latter would be inconsistent.

Similarly, when we call "score" we use the fitted models_nuisance to predict. But which fitted models_nusiance? just the last mc iter? or all the mc_iters?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Right now, it is fitting only a single model multiple times and retaining information from the latest run (because users might expect a single model per fold, not one per fold*mc_iter combination), but I can change that.

Perhaps best would be to have a model wrapper that wraps all of the individual ones but provides a single prediction, similar to bootstrap?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actually I take it back: _crossfit is cloning the model, so we don't call fit repeatedly on the same instance.

@vsyrgkanis
Copy link
Collaborator

@kbattocchi to address #350, we might want to implement public access methods to the cached nuissances, e.g. est.get_nuisances at the ortho_learner level and est.get_residuals at the rlearner level.

Currently this can be done by:

Yres = est._cached_values.nuisances[0]
Tres = est._cached_values.nuisances[0]

which is cumbersome.

return self._discrete_treatment

@discrete_treatment.setter
def discrete_treatment(self, discrete_treatment):
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I found a related bug, not sure if related to the transformer. So the following works:

est = LinearDML(model_y=RandomForestRegressor(),
                model_t=RandomForestClassifier(min_samples_leaf=10),
                discrete_treatment=False,
                linear_first_stages=False,
                n_splits=6)
est.fit(Y, T, X=X, W=W)
te_pred = est.effect(X_test)
te_pred[0]

But the following raises an error, while from the current semantics it's supposed to work fine:

est = LinearDML(model_y=RandomForestRegressor(),
                model_t=RandomForestClassifier(min_samples_leaf=10),
                discrete_treatment=True,
                linear_first_stages=False,
                n_splits=6)
est.fit(Y, T, X=X, W=W)
est.discrete_treatment = False
est.fit(Y, T, X=X, W=W)
te_pred[0]

Which means that there is some state of the object that is maintained but should have been reset when resetting the discrete_treatment.

@kbattocchi kbattocchi force-pushed the kebatt/refit branch 3 times, most recently from 66035e3 to c3a5652 Compare December 23, 2020 16:04
@kbattocchi kbattocchi force-pushed the kebatt/refit branch 2 times, most recently from f09f754 to de96b56 Compare December 23, 2020 19:15
@vsyrgkanis
Copy link
Collaborator

This code fails:

from econml.sklearn_extensions.linear_model import StatsModelsLinearRegression, StatsModelsRLM, DebiasedLasso
from sklearn.linear_model import ElasticNetCV, LassoCV

est = DML(model_y=RandomForestRegressor(min_samples_leaf=20, max_depth=3),
          model_t=RandomForestRegressor(min_samples_leaf=20, max_depth=3),
          model_final=StatsModelsLinearRegression(fit_intercept=False),
          random_state=123)
est.fit(Y_train, T_train, X=X_train, W=W_train, cache_values=True)
print(list(zip(est.coef_, *est.coef__interval())))
est.model_final = StatsModelsRLM(fit_intercept=True)
est.refit()
print(list(zip(est.coef_, *est.coef__interval())))
[(6.45655731600806, 6.050960723017014, 6.8621539089991055)]
4.274931642986886
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
<ipython-input-6-0a8e6a227739> in <module>
     13 est.model_final = StatsModelsRLM(fit_intercept=True)
     14 est.refit()
---> 15 print(list(zip(est.coef_, *est.coef__interval())))
     16 print(est.score_)
     17 

c:\users\vasy\documents\alicedev\econml\econml\cate_estimator.py in call(self, *args, **kwargs)
    205                 return getattr(self._inference, name)(*args, **kwargs)
    206             else:
--> 207                 raise AttributeError("Can't call '%s' because 'inference' is None" % name)
    208         return call
    209 

AttributeError: Can't call 'coef__interval' because 'inference' is None

econml/dml.py Show resolved Hide resolved
econml/dml.py Show resolved Hide resolved
econml/dml.py Show resolved Hide resolved
econml/dml.py Show resolved Hide resolved
# _BaseDML's final model setter reuses its old featurizer
# so we need to pass _BaseDML, not DML to the super() call
# super(_BasaeDML).model_final = _FinalWrapper(...)
super(_BaseDML, _BaseDML).model_final.__set__(self, _FinalWrapper(
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is a bug. There is not model_final

@@ -252,7 +270,7 @@ def model_cate(self):
An instance of the model_final object that was fitted after calling fit which corresponds
to the constant marginal CATE model.
"""
return super().model_final._model
return self._model_final
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This should simply:

self.model_final

without the underscore

@kbattocchi
Copy link
Collaborator Author

Shelving in favor of #360

@kbattocchi kbattocchi closed this Jan 6, 2021
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

4 participants