Can we meaningfully ascribe values to these divergent integrals? My gut feeling is that
$\int_0^\infty (1-\frac1{x^2})dx=0$
$\int_0^\infty (x-\frac2{x^3})dx=0$
$\int_0^\infty (x^2-\frac6{x^4})dx=0,$
etc, and in general,
$\int_0^\infty (x^k-(k+1)!x^{-(k+2)})dx=0,$
is the natural way to regularize these divergent integrals.
The approach is based on Laplace transform (it allows to transform a divergent integral of a function around a pole into a divergent integral of an unbounded function at infinity).
I wonder, whether this approach is known and used anywhere? Does it possess any properties that would justify its "naturalness"?
 A: Depends on how you choose to define it :)
A natural way is to define $\int_0^{\infty} = \lim_{a\rightarrow0} \int_a^{1/a}$, in which case:
\begin{eqnarray*}
\int_0^{\infty}\left(1-\frac1{x^2}\right)dx &=& \lim_{a\rightarrow0} \int_a^{1/a} \left(1-\frac1{x^2}\right)dx \\
&=& \lim_{a\rightarrow0} \left(x+\frac1x\right|_a^{1/a} \\
&=& \lim_{a\rightarrow0} 0 \\
&=& 0
\end{eqnarray*}
If you instead for some reason wanted $\int_0^{\infty} = \lim_{a\rightarrow0} \int_a^{1+1/a}$, this is perfectly valid for convergent integrals and will get the same answer, but for your integral you would now get the value of 1.  Maybe this definition does not look as natural as the original one, but one can easily do natural operations and manipulations on the original expression and end up with a different answer like that.
This is very different from the study of divergent series (see Hardy's book) where you can often show that there's a canonical value you can give to a divergent series (i.e. the now famous $1+1/1+1/3+1/14+\cdots=-1/12$)
A: Warning: The following contains formal manipulations that ignore convergence.
Proposition:
\begin{align}
\mathrm{regularized} \int_0^\infty \mathrm{d}x^s = 0
\end{align}
for all $s$ such that $\Re(s) \neq 0$.

“Proof” 1:
\begin{align}
\int_0^\infty \mathrm{d}x^s
&= \int_0^1 \mathrm{d}x^s + \int_1^\infty \mathrm{d}x^s \\
&= 1 - 1 \\
&= 0
\end{align}
where we simultaneously assumed $\Re(s) > 0$ and $\Re(s) < 0$ for the first and second integrals, respectively.

“Proof” 2:
Let
\begin{align}
f_\pm(\eta) = \mathrm{e}^{-\eta} \left( 1 \pm \frac{\eta}s \right)
\end{align}
Case $\Re(s) > 0$:
\begin{align}
\int_0^\infty x^{s-1} \mathrm{d}x &= \int_0^\infty \lim_{\varepsilon \downarrow 0} x^{s-1} f_-(\varepsilon x) \mathrm{d}x \\
&\stackrel{!}{=} \lim_{\varepsilon \downarrow 0} \int_0^\infty x^{s-1} f_-(\varepsilon x) \mathrm{d}x \\
&= \lim_{\varepsilon \downarrow 0} 0 & \Re(s) > 0, \Re(\varepsilon) > 0 \\
&= 0
\end{align}
Case $\Re(s) < 0$:
\begin{align}
\int_0^\infty x^{s-1} \mathrm{d}x
&= \int_0^\infty \lim_{\varepsilon \downarrow 0} x^{s-1} f_+(\varepsilon x^{-1}) \mathrm{d}x \\
&\stackrel{!}{=} \lim_{\varepsilon \downarrow 0} \int_0^\infty x^{s-1} f_+(\varepsilon x^{-1}) \mathrm{d}x \\
&= \lim_{\varepsilon \downarrow 0} 0 & \Re(s) < 0, \Re(\varepsilon) > 0 \\
&= 0
\end{align}
See here for more context about this regulator-based approach.

Applying linearity of integration yields your equalities.
