3
$\begingroup$

Some background before I state the questions:

I have a $d$-dimensional random vector $X=(X_1,\ldots,X_n)$ and a function $f:\mathbb{R}^d\rightarrow\mathbb{R}$. Ultimately my goal is to understand $f$ and I will explain what I mean by that later. The problem is created by

  1. $d$ being pretty large (e.g. $d>10$)
  2. $f$ being a very complex computer simulation, i.e. no closed formulas exist and the computation of $f$ is expensive in terms of time and resources

Since samples from $X$ are comparatively cheap, a standard approach is to approximate $f$ by a parametric function $g$ which is easier to handle. In practice, parameters of $g$ are estimated using the values of $f$ on a suitable sample from $X$ using least squares.

Note that this is different from regression, since there is NO additive stochastic error term, i.e. $f$ and $g$ are both deterministic functions of $X$. On the other hand it seems to differ from (plain vanilla) function approximation, since I have an underlying probability space and need to do some sort of statistical estimation.

So my first question is: It is not regression, it is not approximation, so what is it?

Some more specific questions:

It is tempting to use the regression framework anyway, since this is a convenient way to produce the OLS estimates for the parameters. It seems that all inference and diagnostics derived from properties of the error distribution in regression will be meaningless (e.g. F-tests for the significance of coefficients, adjusted R squared and so on). Is there anything useful, beyond the parameter estimate itself, which survives from the regression framework?

What type of estimation is possible in this context? For example, I might assume that $f$ is a sum of univariate functions (i.e. only main effects). Are there tests to see whether this is the case or whether I need interactions? Another typical question: Are all $X_i$ equally important and what are the least/most important inputs? Those questions have answers in the context of regression. What is possible here?

Of course the answers will depend on assumptions for $f$. What would be useful ones? Actually smoothness properties are not entirely clear but would for example piecewise linear be helpful?

$\endgroup$
  • $\begingroup$ How expensive is each function evaluation? That is, how many "observations" (function evaluation) can you take? About 10 (then the problem becomes almost impossible), about 100? $\endgroup$ – kjetil b halvorsen Oct 13 '14 at 14:59
  • $\begingroup$ A few thousand should be possible. But this sounds better than it is since I am mostly interested in tail behaviour (99% say). $\endgroup$ – g g Oct 13 '14 at 17:19
3
$\begingroup$

$\newcommand{\R}{\mathbb{R}}$ You can use some kind of regression here! You have a (complicated, costly) model $y=f(x)$ that you want to approximate, where $x \in \R^d$. You can define some kind of regression model (linear or nonlinear, I will use a linear model in the discussion) to approximate $f$.

To take into account interactions and maybe other nonlinear transformation of variables, define an observation vector $x^*$ which is some (known) function of $x$, which have as components the components of $x$, products of components of $x$ (interaction terms) and so on. That gives as model $y=\beta^T x^*$ for some parameter vector $\beta$ to estimate. The error in this regression model is $f(x)-\beta^T x^*$.

What is different from the usual regression model is that this error function in principle is a known, deterministic function, (only to costly to calculate). So this does not make sense as a statistical model in the usual frequentist sense, since the errors will be identical on repeated application, so do not have a distribution in the frequentist sense. But we can make sense of this model in bayesian terms, so the error term can have a distribution which represents our lack of knowledge. But typically, errors for nearby points (in the space of $x$-values) will be similar, so we need a model for dependent errors. One possibility is to assume the errors form a stationary Gaussian random field. Regression under such a model is known as kriging, and heavily used in spatial statistics. Kriging outside spatial statistics is also known as gaussian process models.

The use of statistical models as above to represent the output of complicated computer models, is known in some circles as metamodelling. This could be seen as a kind of bayesian modeling, the metamodel working as a kind of prior.

Some references is given below:

http://amstat.tandfonline.com/doi/abs/10.1080/01621459.1991.10475138#.VD4ySvl_t8E

http://mechanicaldesign.asmedigitalcollection.asme.org/article.aspx?articleid=1449318

$\endgroup$
  • $\begingroup$ I just noticed I can't vote up (not enough reputation). Sorry for that because the answer and references seem to be very helpful. $\endgroup$ – g g Oct 15 '14 at 9:53
  • $\begingroup$ You can come back and vote when you get more reputationn points... $\endgroup$ – kjetil b halvorsen Oct 15 '14 at 9:59
  • $\begingroup$ you might want to search for gaussian process optimization eg robots.ox.ac.uk/~parg/pubs/OsborneGarnettRobertsGPGO.pdf - assuming the aim of using this approximation is to drive optimisation... the ideas here being that you are explicitly modelling which parts of the function you have more certainty about and which may be worth exploring $\endgroup$ – seanv507 Oct 15 '14 at 10:06
1
$\begingroup$

Note that this is different from regression, since there is NO additive error term, i.e. f and g are both deterministic functions of X. On the other hand it seems to differ from (plain vanilla) function approximation, since I have an underlying probability space and need to do some sort of statistical estimation.

Actually there is always an additive error term. By definition it's the difference between $f$ and $g$ at some value of $X$. Perhaps you mean that it's not distributed according to some convenient function, or not constant across $X$, or something like that. All those may be true, but none of them stop you using a regression framework - that is, thinking of $g$ as a model of the expected value of $f$, i.e. a model of $E[f(X)]$. This has a lot of optimality properties even without assuming a convenient distribution for the errors, as any econometrics textbook will expound at length. Which regression framework you use is a further, more interesting question.

For your purposes the linear additive approach of basic textbooks is probably not going to do so well unless your $X$'s are very carefully chosen. The non-linear extensions of this approach e.g. Gaussian Processes, splines, GAMs etc. assume that $E[f(X)]$ changes smoothly (if not necessarily linearly) with changes in X. This might work for you, but I'd also guess not.

A potentially interesting alternative would be to use some kind of decision tree approach e.g. CART, which recursively partitions X on the assumption that discontinuities in $E[f(X)]$ occur over boundaries in elements of X.

The limiting case of both of these approaches is a saturated dummy variable regression model that is, on the one hand linear again, but on the other, requires more data than you have got to fit.

The best advice is to try a bunch and maybe average their predictions. Performance will depend on how well the bias of the methods you use e.g. X-aligned discontinuities, vs. smoothness, fit the actual $X\rightarrow E[f(X)]$ function. But there's nothing in principle problematic in thinking of this as a regression problem.

$\endgroup$
  • $\begingroup$ Maybe I was not precise in my terminology. I should have said additive stochastic(!) error term. The difference between $f$ and $g$ at any value of $X$ is NOT a random variable in my setup! So it is VERY conveniently distributed: It is just a constant number hence $E[f(X)\mid X=x]=f(x)$ and nothing is gained. $\endgroup$ – g g Oct 15 '14 at 17:51
  • $\begingroup$ I think it is rather that I was being imprecise: unless X spans all possible system inputs (which it might) and your $g$ is exactly $f$ then you will have residual uncertainty about $f$ at $x$ on the basis of a fitted model value $g(x)$, and it will be summarised by the difference between them. This difference will be because of the imperfection of $g$ as a model of $f$. So that's your residual. And it has a distribution too, even when $f(X)$ is a deterministic function of $X$ (unless the model is perfect). I may do a little answer adjusting to reflect that... $\endgroup$ – conjugateprior Oct 16 '14 at 8:00
  • $\begingroup$ I think the distinction between the "standard regression approach" and the approach suggested by kjetil above is important. Do you agree that $f(X)=g(X) + \epsilon$, where $\epsilon$ is a random variable is a "standard regression approach"? Hence conditional on a single $X=x$ the response $f(x)$ needs to be a random variable. Do you agree that this is not so in my case? I found this distinction hard to see at first. This is why I have slight trouble using the term "regression" for both and insist a bit on clarity. $\endgroup$ – g g Oct 16 '14 at 9:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.