Predicting values in a given range [0, 1] (not probability) I need to predict the impact of a set of node failures in a network, based on 2 features: the fraction of failed nodes and a measure of their network centrality.
Failure of less important nodes will have low to no repercussions, while failure of important nodes can bring down the whole network.
The predicted output should be a number in the range [0, 1], where 0 is "no damage", and "1" is "failure of all nodes". 0 is inclusive since I use it to cover the case "no damage" when the fraction of failed nodes is 0.
Here's an example mapping 2 features (x and y axes) to the actual simulation result (z axis). I plotted the scatter and the projection of the points. It looks like something that can be learned, right?


I'm using multivariate linear regression in scikit-learn to train a predictor. I take my training set, pick the 2 features I want, use them to make a 4th degree polynomial, which I then standardize and use to train my model.
The problem is, for some combinations of my 2 feature values, the predictor outputs negative values or values greater than 1. I can post a picture if needed.
I've tried changing the degree of the polynomial, using regularization (lasso, ridge and ElasticNet) and applying normalization. Nothing completely fixes the problem, especially for predictions >1.
Is linear regression the wrong tool for predicting values in a range?
Any suggestion for a drop-in replacement?
I'm tempted to try logistic regression, but it does not sound like the right tool, since it's for classification, and it outputs probabilities, while my values have a different meaning (not a "chance" but a degree of damage).
 A: You are facing extremely different responses to very similar inputs. Regression will not give you meaningful answers.
Classification approach.
Consider, alternatively, categorical responses such as 'low damage' and 'high damage'. You should be looking for a function $f$ such that $f(p, c)$ is the probability that you will have 'low damage' when a fraction $p$ of the nodes have failed with loss of centrality $c$.
You can construct $f$ in the following way. Let $n \in \mathbb{N}$, let $a_{i,j}$ be the number of points in the region $R_{i,j} = [\frac{i}{n+1}, \frac{i+1}{n+1}]\times [\frac{j}{n+1}, \frac{j+1}{n+1}]$ corresponding to 'low damage' and let $N_{i,j}$ be the total number of points in this region. Thus
$$
f_{i,j} = \frac{a_{i,j}}{N_{i,j}}
$$
is the fraction of points in $R_{i,j}$ corresponding to 'low damage'. You can then let
$$
f(p,c) = f_{i,j}\quad  \text{ if } (p,c) \in R_{i,j}.
$$
If you want something smooth, you can use
$$
f(p,c) = \sum_{i,j} f_{i,j} B_{i,n}(p) B_{j,n}(c),
$$
where $B_{i,n}(p) = {n \choose i} x^{i}(1-x)^{n-i}$ are the Bernstein polynomials.
Improving the model.
You could put a prior on the $f_{i,j}$ and on $n$. This would give you a reasonable bayesian model, and fix the small problem that arises when $N_{i,j} = 0$.
Visualizing the results.
The set of points $(p,c)$ such that $f(p,c) = \frac{1}{2}$ will give you a classification boundary that is an algebraic variety. You also have critical regions $C_\varepsilon = \{(p,c) \,|\, f(p,c) < \varepsilon\}$ where 'high damage' is very likely. 
You can calculate the marginal probability of high or low damages for a given fraction of failed nodes and plot this together with the probability surface.
