Repository: pydata/xarray

Issue:
Bug in WarpedVRT support of open_rasterio()
#### Code Sample, a copy-pastable example if possible

Using the following data:

https://gitlab.orfeo-toolbox.org/orfeotoolbox/otb/blob/develop/Data/Input/QB_Toulouse_Ortho_XS.tif

```python
import rasterio as rio
from rasterio.crs import CRS
from rasterio.warp import calculate_default_transform,aligned_target
from rasterio.enums import Resampling
from rasterio.vrt import WarpedVRT

import xarray as xr

path = 'QB_Toulouse_Ortho_XS.tif'

# Read input metadata with rasterio                                                                                                                                                                                                           
with rio.open(path) as src:
    print('Input file CRS is {}'.format(src.crs))
    print('Input file shape is {}'.format(src.shape))
    print('Input file transform is {}'.format(src.transform))
    # Create a different CRS                                                                                                                                                                                                                  
    dst_crs = CRS.from_epsg(2154)
    # Compute a transform that will resample to dst_crs and change resolution to dst_crs                                                                                                                                                      
    transform, width, height = calculate_default_transform(
        src.crs, dst_crs, src.width, src.height,resolution=20, *src.bounds)
    # Fill vrt options as shown in https://rasterio.readthedocs.io/en/stable/topics/virtual-warping.html                                                                                                                                      
    vrt_options = {
    'resampling': Resampling.cubic,
    'crs': dst_crs,
    'transform': transform,
    'height': height,
    'width': width
    }
    # Create a WarpedVRT using the vrt_options                                                                                                                                                                                                
    with WarpedVRT(src,**vrt_options) as vrt:
        print('VRT shape is {}'.format(vrt.shape))
        # Open VRT with xarray                                                                                                                                                                                                                
        ds = xr.open_rasterio(vrt)
        # Shape does not match vrt shape:                                                                                                                                                                                                     
        print(ds)
```
Output:
```
$ python test_rio_vrt.py
Input file CRS is EPSG:32631
Input file shape is (500, 500)
Input file transform is | 0.60, 0.00, 374149.98|
| 0.00,-0.60, 4829183.99|
| 0.00, 0.00, 1.00|
VRT shape is (16, 16)
<xarray.DataArray (band: 4, y: 500, x: 500)>
[1000000 values with dtype=int16]
Coordinates:
  * band     (band) int64 1 2 3 4
  * y        (y) float64 6.28e+06 6.28e+06 6.28e+06 ... 6.279e+06 6.279e+06
  * x        (x) float64 5.741e+05 5.741e+05 5.741e+05 ... 5.744e+05 5.744e+05
Attributes:
    transform:   (0.6003151072155879, 0.0, 574068.2261249251, 0.0, -0.6003151...
    crs:         EPSG:2154
    res:         (0.6003151072155879, 0.6003151072155879)
    is_tiled:    0
    nodatavals:  (nan, nan, nan, nan)
```

#### Problem description

In the above example, `xarray.open_rasterio()` is asked to read a `WarpedVRT` created with `rasterio`. This `WarpedVRT` has custom `transform`, `width` and `height`, which is a very common use case of `WarpedVRT` where you upsample or downsample your data on the fly during read. Las, `open_rasterio()` ignores those custom attributes, which result in a wrong `DataArray` shape:
```
VRT shape is (16, 16)
<xarray.DataArray (band: 4, y: 500, x: 500)>
```

#### Expected Output
Correct output should be:

```
VRT shape is (16, 16)
<xarray.DataArray (band: 4, y: 16, x: 16)>
```

Which is fairly easy to obtain by modifying the following lines in `open_rasterio()`:
https://github.com/pydata/xarray/blob/0c73a380745c4792ab440eb020f78f203897abe5/xarray/backends/rasterio_.py#L222

With the following:
```python
      vrt_params = dict(crs=vrt.crs.to_string(),
                          resampling=vrt.resampling,
                          src_nodata=vrt.src_nodata,
                          dst_nodata=vrt.dst_nodata,
                          tolerance=vrt.tolerance,
                          # Edit                                                                                                                                                                                                              
                          transform=vrt.transform,
                          width=vrt.width,
                          height=vrt.height,
                          # End edit                                                                                                                                                                                                          
                          warp_extras=vrt.warp_extras)
```
I can provide this patch in a pull request if needed.

#### Output of ``xr.show_versions()``

<details>
>>> xr.show_versions()

INSTALLED VERSIONS
------------------
commit: None
python: 3.6.7 | packaged by conda-forge | (default, Nov 21 2018, 03:09:43) 
[GCC 7.3.0]
python-bits: 64
OS: Linux
OS-release: 3.10.0-862.2.3.el7.x86_64
machine: x86_64
processor: x86_64
byteorder: little
LC_ALL: None
LANG: fr_FR.utf8
LOCALE: fr_FR.UTF-8
libhdf5: 1.10.4
libnetcdf: 4.6.2

xarray: 0.11.3
pandas: 0.24.1
numpy: 1.16.1
scipy: 1.2.0
netCDF4: 1.4.2
pydap: None
h5netcdf: None
h5py: None
Nio: None
zarr: None
cftime: 1.0.3.4
PseudonetCDF: None
rasterio: 1.0.17
cfgrib: None
iris: None
bottleneck: None
cyordereddict: None
dask: 1.1.1
distributed: 1.25.3
matplotlib: 3.0.2
cartopy: 0.17.0
seaborn: 0.9.0
setuptools: 40.7.3
pip: 19.0.1
conda: None
pytest: 4.2.0
IPython: 7.1.1
sphinx: None

</details>


Discussion:
fmaussion: >  Should I add tests for this ?

If there is a way, it would be great yes. You can build upon (or adapt) this example here:

https://github.com/pydata/xarray/blob/225868d232219440b188956531d5764ff4cd1b53/xarray/tests/test_backends.py#L3327-L3351

cc @scottyhq who wrote this code, as I don't think any of the xarray team uses this format.

Please also add a line to [whats-new.rst](https://github.com/pydata/xarray/blob/master/doc/whats-new.rst)

pep8speaks: Hello @jmichel-otb! Thanks for updating this PR. We checked the lines you've touched for [PEP 8](https://www.python.org/dev/peps/pep-0008) issues, and found:



* In the file [`xarray/tests/test_backends.py`](https://github.com/pydata/xarray/blob/f8de48c7a606e4b0209388bac4ec01e343dfc3d2/xarray/tests/test_backends.py):

> [Line 3359:80](https://github.com/pydata/xarray/blob/f8de48c7a606e4b0209388bac4ec01e343dfc3d2/xarray/tests/test_backends.py#L3359): [E501](https://duckduckgo.com/?q=pep8%20E501) line too long (85 > 79 characters)
> [Line 3361:80](https://github.com/pydata/xarray/blob/f8de48c7a606e4b0209388bac4ec01e343dfc3d2/xarray/tests/test_backends.py#L3361): [E501](https://duckduckgo.com/?q=pep8%20E501) line too long (85 > 79 characters)
> [Line 3362:80](https://github.com/pydata/xarray/blob/f8de48c7a606e4b0209388bac4ec01e343dfc3d2/xarray/tests/test_backends.py#L3362): [E501](https://duckduckgo.com/?q=pep8%20E501) line too long (89 > 79 characters)
> [Line 3363:80](https://github.com/pydata/xarray/blob/f8de48c7a606e4b0209388bac4ec01e343dfc3d2/xarray/tests/test_backends.py#L3363): [E501](https://duckduckgo.com/?q=pep8%20E501) line too long (105 > 79 characters)


jmichel-otb: @fmaussion done, let's see what CI has to say about my patches ;)

I remember reading a thread somewhere on xarray github repo discussing whether xarray should include the rasterio backend or not.

I understand that bridges between two libraries are always hard to maintain, because you need to know both products (we actually have the same kind of problem with OTB and QGis), but from a user standpoint, they need to exist somewhere. I would probably never have turned to xarray if someone with the required knowledge had not implemented the rasterio backend.

Then of course the user community should take care of maintaining those backends (this is what I am doing right now).

Bridging xarray with rasterio opens xarray to the remote sensing imagery community. And behind rasterio there is gdal, which is an awesome library with so many great capabilities (like this on-the-fly reprojection during reading I mentioned).

scottyhq: @jmichel-otb completely agree on your comments and think that the remote sensing community gains a lot by making these packages work well together (gdal, rasterio, xarray, dask, etc). As more people start using the rasterio backend with various use cases hopefully there will be more contributions such as this. Thanks!

Patch:
diff --git a/.travis.yml b/.travis.yml
index ea9ee7ad..212ddb77 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -65,7 +65,7 @@ script:
     elif [[ "$CONDA_ENV" == "py36-hypothesis" ]]; then
       pytest properties ;
     else
-      py.test xarray --cov=xarray --cov-config ci/.coveragerc --cov-report term-missing --verbose $EXTRA_FLAGS;
+      py.test xarray --cov=xarray --cov-config ci/.coveragerc --cov-report term-missing $EXTRA_FLAGS;
     fi
 
 after_success:
diff --git a/README.rst b/README.rst
index 6dbf7745..83382f87 100644
--- a/README.rst
+++ b/README.rst
@@ -8,9 +8,9 @@ xarray: N-D labeled arrays and datasets
 .. image:: https://coveralls.io/repos/pydata/xarray/badge.svg
    :target: https://coveralls.io/r/pydata/xarray
 .. image:: https://readthedocs.org/projects/xray/badge/?version=latest
-   :target: http://xarray.pydata.org/
-.. image:: http://img.shields.io/badge/benchmarked%20by-asv-green.svg?style=flat
-  :target: http://pandas.pydata.org/speed/xarray/
+   :target: https://xarray.pydata.org/
+.. image:: https://img.shields.io/badge/benchmarked%20by-asv-green.svg?style=flat
+  :target: https://pandas.pydata.org/speed/xarray/
 .. image:: https://img.shields.io/pypi/v/xarray.svg
    :target: https://pypi.python.org/pypi/xarray/
 
@@ -30,10 +30,10 @@ It is particularly tailored to working with netCDF_ files, which were the
 source of xarray's data model, and integrates tightly with dask_ for parallel
 computing.
 
-.. _NumPy: http://www.numpy.org
-.. _pandas: http://pandas.pydata.org
-.. _dask: http://dask.org
-.. _netCDF: http://www.unidata.ucar.edu/software/netcdf
+.. _NumPy: https://www.numpy.org
+.. _pandas: https://pandas.pydata.org
+.. _dask: https://dask.org
+.. _netCDF: https://www.unidata.ucar.edu/software/netcdf
 
 Why xarray?
 -----------
@@ -66,12 +66,12 @@ powerful and concise interface. For example:
 Documentation
 -------------
 
-Learn more about xarray in its official documentation at http://xarray.pydata.org/
+Learn more about xarray in its official documentation at https://xarray.pydata.org/
 
 Contributing
 ------------
 
-You can find information about contributing to xarray at our `Contributing page <http://xarray.pydata.org/en/latest/contributing.html#>`_.
+You can find information about contributing to xarray at our `Contributing page <https://xarray.pydata.org/en/latest/contributing.html#>`_.
 
 Get in touch
 ------------
@@ -81,9 +81,9 @@ Get in touch
 - For less well defined questions or ideas, or to announce other projects of
   interest to xarray users, use the `mailing list`_.
 
-.. _StackOverFlow: http://stackoverflow.com/questions/tagged/python-xarray
+.. _StackOverFlow: https://stackoverflow.com/questions/tagged/python-xarray
 .. _mailing list: https://groups.google.com/forum/#!forum/xarray
-.. _on GitHub: http://github.com/pydata/xarray
+.. _on GitHub: https://github.com/pydata/xarray
 
 NumFOCUS
 --------
@@ -120,7 +120,7 @@ Licensed under the Apache License, Version 2.0 (the "License");
 you may not use this file except in compliance with the License.
 You may obtain a copy of the License at
 
-  http://www.apache.org/licenses/LICENSE-2.0
+  https://www.apache.org/licenses/LICENSE-2.0
 
 Unless required by applicable law or agreed to in writing, software
 distributed under the License is distributed on an "AS IS" BASIS,
diff --git a/doc/whats-new.rst b/doc/whats-new.rst
index 9c88445b..692082b1 100644
--- a/doc/whats-new.rst
+++ b/doc/whats-new.rst
@@ -13,26 +13,40 @@ What's New
     import xarray as xr
     np.random.seed(123456)
 
-.. _whats-new.0.12.1:
+.. _whats-new.0.12.2:
 
-v0.12.1 (unreleased)
+v0.12.2 (unreleased)
 --------------------
 
 Enhancements
 ~~~~~~~~~~~~
 
+Bug fixes
+~~~~~~~~~
+
+.. _whats-new.0.12.1:
+
+v0.12.1 (4 April 2019)
+----------------------
+
+Enhancements
+~~~~~~~~~~~~
+
 - Allow ``expand_dims`` method to support inserting/broadcasting dimensions
   with size > 1. (:issue:`2710`)
   By `Martin Pletcher <https://github.com/pletchm>`_.
 
-
 Bug fixes
 ~~~~~~~~~
 
 - Dataset.copy(deep=True) now creates a deep copy of the attrs (:issue:`2835`).
   By `Andras Gefferth <https://github.com/kefirbandi>`_.
-- ``swap_dims`` would create incorrect ``indexes`` (:issue:`2842`).
+- Fix incorrect ``indexes`` resulting from various ``Dataset`` operations
+  (e.g., ``swap_dims``, ``isel``, ``reindex``, ``[]``) (:issue:`2842`,
+  :issue:`2856`).
   By `Stephan Hoyer <https://github.com/shoyer>`_.
+- open_rasterio() now supports rasterio.vrt.WarpedVRT with custom transform, width and height (:issue:`2864`).
+  By `Julien Michel <https://github.com/jmichel-otb>`_.
 
 .. _whats-new.0.12.0:
 
diff --git a/xarray/backends/rasterio_.py b/xarray/backends/rasterio_.py
index b629eb51..2f00702f 100644
--- a/xarray/backends/rasterio_.py
+++ b/xarray/backends/rasterio_.py
@@ -224,6 +224,9 @@ def open_rasterio(filename, parse_coordinates=None, chunks=None, cache=None,
                           src_nodata=vrt.src_nodata,
                           dst_nodata=vrt.dst_nodata,
                           tolerance=vrt.tolerance,
+                          transform=vrt.transform,
+                          width=vrt.width,
+                          height=vrt.height,
                           warp_extras=vrt.warp_extras)
 
     if lock is None:
diff --git a/xarray/coding/cftime_offsets.py b/xarray/coding/cftime_offsets.py
index d724554b..2ee38a20 100644
--- a/xarray/coding/cftime_offsets.py
+++ b/xarray/coding/cftime_offsets.py
@@ -79,7 +79,7 @@ def get_date_type(calendar):
 
 class BaseCFTimeOffset(object):
     _freq = None  # type: ClassVar[str]
-    _day_option = None
+    _day_option = None  # type: ClassVar[str]
 
     def __init__(self, n=1):
         if not isinstance(n, int):
diff --git a/xarray/core/alignment.py b/xarray/core/alignment.py
index af08eef2..642be735 100644
--- a/xarray/core/alignment.py
+++ b/xarray/core/alignment.py
@@ -315,36 +315,51 @@ def reindex_variables(
     """
     from .dataarray import DataArray
 
+    # create variables for the new dataset
+    reindexed = OrderedDict()  # type: OrderedDict[Any, Variable]
+
     # build up indexers for assignment along each dimension
     int_indexers = {}
-    targets = OrderedDict()  # type: OrderedDict[Any, pd.Index]
+    new_indexes = OrderedDict(indexes)
     masked_dims = set()
     unchanged_dims = set()
 
-    # size of reindexed dimensions
-    new_sizes = {}
+    for dim, indexer in indexers.items():
+        if isinstance(indexer, DataArray) and indexer.dims != (dim,):
+            warnings.warn(
+                "Indexer has dimensions {0:s} that are different "
+                "from that to be indexed along {1:s}. "
+                "This will behave differently in the future.".format(
+                    str(indexer.dims), dim),
+                FutureWarning, stacklevel=3)
+
+        target = new_indexes[dim] = utils.safe_cast_to_index(indexers[dim])
+
+        if dim in indexes:
+            index = indexes[dim]
 
-    for name, index in indexes.items():
-        if name in indexers:
             if not index.is_unique:
                 raise ValueError(
                     'cannot reindex or align along dimension %r because the '
-                    'index has duplicate values' % name)
-
-            target = utils.safe_cast_to_index(indexers[name])
-            new_sizes[name] = len(target)
+                    'index has duplicate values' % dim)
 
             int_indexer = get_indexer_nd(index, target, method, tolerance)
 
             # We uses negative values from get_indexer_nd to signify
             # values that are missing in the index.
             if (int_indexer < 0).any():
-                masked_dims.add(name)
+                masked_dims.add(dim)
             elif np.array_equal(int_indexer, np.arange(len(index))):
-                unchanged_dims.add(name)
+                unchanged_dims.add(dim)
 
-            int_indexers[name] = int_indexer
-            targets[name] = target
+            int_indexers[dim] = int_indexer
+
+        if dim in variables:
+            var = variables[dim]
+            args = (var.attrs, var.encoding)  # type: tuple
+        else:
+            args = ()
+        reindexed[dim] = IndexVariable((dim,), target, *args)
 
     for dim in sizes:
         if dim not in indexes and dim in indexers:
@@ -356,25 +371,6 @@ def reindex_variables(
                     'index because its size %r is different from the size of '
                     'the new index %r' % (dim, existing_size, new_size))
 
-    # create variables for the new dataset
-    reindexed = OrderedDict()  # type: OrderedDict[Any, Variable]
-
-    for dim, indexer in indexers.items():
-        if isinstance(indexer, DataArray) and indexer.dims != (dim,):
-            warnings.warn(
-                "Indexer has dimensions {0:s} that are different "
-                "from that to be indexed along {1:s}. "
-                "This will behave differently in the future.".format(
-                    str(indexer.dims), dim),
-                FutureWarning, stacklevel=3)
-
-        if dim in variables:
-            var = variables[dim]
-            args = (var.attrs, var.encoding)  # type: tuple
-        else:
-            args = ()
-        reindexed[dim] = IndexVariable((dim,), indexers[dim], *args)
-
     for name, var in variables.items():
         if name not in indexers:
             key = tuple(slice(None)
@@ -395,9 +391,6 @@ def reindex_variables(
 
             reindexed[name] = new_var
 
-    new_indexes = OrderedDict(indexes)
-    new_indexes.update(targets)
-
     return reindexed, new_indexes
 
 
diff --git a/xarray/core/common.py b/xarray/core/common.py
index 6ec07156..b518e843 100644
--- a/xarray/core/common.py
+++ b/xarray/core/common.py
@@ -1,6 +1,8 @@
 from collections import OrderedDict
 from contextlib import suppress
 from textwrap import dedent
+from typing import (Any, Callable, Hashable, Iterable, Iterator, List, Mapping,
+                    MutableMapping, Optional, Tuple, TypeVar, Union)
 
 import numpy as np
 import pandas as pd
@@ -11,13 +13,18 @@ from .options import _get_keep_attrs
 from .pycompat import dask_array_type
 from .utils import Frozen, ReprObject, SortedKeysDict, either_dict_or_kwargs
 
+
 # Used as a sentinel value to indicate a all dimensions
 ALL_DIMS = ReprObject('<all-dims>')
 
 
-class ImplementsArrayReduce(object):
+T = TypeVar('T')
+
+
+class ImplementsArrayReduce:
     @classmethod
-    def _reduce_method(cls, func, include_skipna, numeric_only):
+    def _reduce_method(cls, func: Callable, include_skipna: bool,
+                       numeric_only: bool):
         if include_skipna:
             def wrapped_func(self, dim=None, axis=None, skipna=None,
                              **kwargs):
@@ -46,9 +53,10 @@ class ImplementsArrayReduce(object):
             and 'axis' arguments can be supplied.""")
 
 
-class ImplementsDatasetReduce(object):
+class ImplementsDatasetReduce:
     @classmethod
-    def _reduce_method(cls, func, include_skipna, numeric_only):
+    def _reduce_method(cls, func: Callable, include_skipna: bool,
+                       numeric_only: bool):
         if include_skipna:
             def wrapped_func(self, dim=None, skipna=None,
                              **kwargs):
@@ -76,46 +84,38 @@ class ImplementsDatasetReduce(object):
 
 
 class AbstractArray(ImplementsArrayReduce):
-    """Shared base class for DataArray and Variable."""
-
-    def __bool__(self):
+    """Shared base class for DataArray and Variable.
+    """
+    def __bool__(self: Any) -> bool:
         return bool(self.values)
 
-    # Python 3 uses __bool__, Python 2 uses __nonzero__
-    __nonzero__ = __bool__
-
-    def __float__(self):
+    def __float__(self: Any) -> float:
         return float(self.values)
 
-    def __int__(self):
+    def __int__(self: Any) -> int:
         return int(self.values)
 
-    def __complex__(self):
+    def __complex__(self: Any) -> complex:
         return complex(self.values)
 
-    def __long__(self):
-        return long(self.values)  # noqa
-
-    def __array__(self, dtype=None):
+    def __array__(self: Any, dtype: Union[str, np.dtype, None] = None
+                  ) -> np.ndarray:
         return np.asarray(self.values, dtype=dtype)
 
-    def __repr__(self):
+    def __repr__(self) -> str:
         return formatting.array_repr(self)
 
-    def _iter(self):
+    def _iter(self: Any) -> Iterator[Any]:
         for n in range(len(self)):
             yield self[n]
 
-    def __iter__(self):
+    def __iter__(self: Any) -> Iterator[Any]:
         if self.ndim == 0:
             raise TypeError('iteration over a 0-d array')
         return self._iter()
 
-    @property
-    def T(self):
-        return self.transpose()
-
-    def get_axis_num(self, dim):
+    def get_axis_num(self, dim: Union[Hashable, Iterable[Hashable]]
+                     ) -> Union[int, Tuple[int, ...]]:
         """Return axis number(s) corresponding to dimension(s) in this array.
 
         Parameters
@@ -128,12 +128,12 @@ class AbstractArray(ImplementsArrayReduce):
         int or tuple of int
             Axis number or numbers corresponding to the given dimensions.
         """
-        if isinstance(dim, str):
-            return self._get_axis_num(dim)
-        else:
+        if isinstance(dim, Iterable) and not isinstance(dim, str):
             return tuple(self._get_axis_num(d) for d in dim)
+        else:
+            return self._get_axis_num(dim)
 
-    def _get_axis_num(self, dim):
+    def _get_axis_num(self: Any, dim: Hashable) -> int:
         try:
             return self.dims.index(dim)
         except ValueError:
@@ -141,7 +141,7 @@ class AbstractArray(ImplementsArrayReduce):
                              (dim, self.dims))
 
     @property
-    def sizes(self):
+    def sizes(self: Any) -> Mapping[Hashable, int]:
         """Ordered mapping from dimension names to lengths.
 
         Immutable.
@@ -153,7 +153,7 @@ class AbstractArray(ImplementsArrayReduce):
         return Frozen(OrderedDict(zip(self.dims, self.shape)))
 
 
-class AttrAccessMixin(object):
+class AttrAccessMixin:
     """Mixin class that allows getting keys with attribute access
     """
     _initialized = False
@@ -168,7 +168,7 @@ class AttrAccessMixin(object):
         """List of places to look-up items for key-autocompletion """
         return []
 
-    def __getattr__(self, name):
+    def __getattr__(self, name: str) -> Any:
         if name != '__setstate__':
             # this avoids an infinite loop when pickle looks for the
             # __setstate__ attribute before the xarray object is initialized
@@ -178,7 +178,7 @@ class AttrAccessMixin(object):
         raise AttributeError("%r object has no attribute %r" %
                              (type(self).__name__, name))
 
-    def __setattr__(self, name, value):
+    def __setattr__(self, name: str, value: Any) -> None:
         if self._initialized:
             try:
                 # Allow setting instance variables if they already exist
@@ -192,7 +192,7 @@ class AttrAccessMixin(object):
                     "assign variables." % (name, type(self).__name__))
         object.__setattr__(self, name, value)
 
-    def __dir__(self):
+    def __dir__(self) -> List[str]:
         """Provide method name lookup and completion. Only provide 'public'
         methods.
         """
@@ -202,7 +202,7 @@ class AttrAccessMixin(object):
                        if isinstance(item, str)]
         return sorted(set(dir(type(self)) + extra_attrs))
 
-    def _ipython_key_completions_(self):
+    def _ipython_key_completions_(self) -> List[str]:
         """Provide method for the key-autocompletions in IPython.
         See http://ipython.readthedocs.io/en/stable/config/integrating.html#tab-completion
         For the details.
@@ -214,49 +214,57 @@ class AttrAccessMixin(object):
         return list(set(item_lists))
 
 
-def get_squeeze_dims(xarray_obj, dim, axis=None):
+def get_squeeze_dims(xarray_obj,
+                     dim: Union[Hashable, Iterable[Hashable], None] = None,
+                     axis: Union[int, Iterable[int], None] = None
+                     ) -> List[Hashable]:
     """Get a list of dimensions to squeeze out.
     """
     if dim is not None and axis is not None:
         raise ValueError('cannot use both parameters `axis` and `dim`')
-
     if dim is None and axis is None:
-        dim = [d for d, s in xarray_obj.sizes.items() if s == 1]
+        return [d for d, s in xarray_obj.sizes.items() if s == 1]
+
+    if isinstance(dim, Iterable) and not isinstance(dim, str):
+        dim = list(dim)
+    elif dim is not None:
+        dim = [dim]
     else:
-        if isinstance(dim, str):
-            dim = [dim]
+        assert axis is not None
         if isinstance(axis, int):
-            axis = (axis, )
-        if isinstance(axis, tuple):
-            for a in axis:
-                if not isinstance(a, int):
-                    raise ValueError(
-                        'parameter `axis` must be int or tuple of int.')
-            alldims = list(xarray_obj.sizes.keys())
-            dim = [alldims[a] for a in axis]
-        if any(xarray_obj.sizes[k] > 1 for k in dim):
-            raise ValueError('cannot select a dimension to squeeze out '
-                             'which has length greater than one')
+            axis = [axis]
+        axis = list(axis)
+        if any(not isinstance(a, int) for a in axis):
+            raise TypeError(
+                'parameter `axis` must be int or iterable of int.')
+        alldims = list(xarray_obj.sizes.keys())
+        dim = [alldims[a] for a in axis]
+
+    if any(xarray_obj.sizes[k] > 1 for k in dim):
+        raise ValueError('cannot select a dimension to squeeze out '
+                         'which has length greater than one')
     return dim
 
 
 class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
     """Shared base class for Dataset and DataArray."""
 
-    def squeeze(self, dim=None, drop=False, axis=None):
+    def squeeze(self, dim: Union[Hashable, Iterable[Hashable], None] = None,
+                drop: bool = False,
+                axis: Union[int, Iterable[int], None] = None):
         """Return a new object with squeezed data.
 
         Parameters
         ----------
-        dim : None or str or tuple of str, optional
+        dim : None or Hashable or iterable of Hashable, optional
             Selects a subset of the length one dimensions. If a dimension is
             selected with length greater than one, an error is raised. If
             None, all length one dimensions are squeezed.
         drop : bool, optional
             If ``drop=True``, drop squeezed coordinates instead of making them
             scalar.
-        axis : int, optional
-            Select the dimension to squeeze. Added for compatibility reasons.
+        axis : None or int or iterable of int, optional
+            Like dim, but positional.
 
         Returns
         -------
@@ -271,7 +279,7 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
         dims = get_squeeze_dims(self, dim, axis)
         return self.isel(drop=drop, **{d: 0 for d in dims})
 
-    def get_index(self, key):
+    def get_index(self, key: Hashable) -> pd.Index:
         """Get an index for a dimension, with fall-back to a default RangeIndex
         """
         if key not in self.dims:
@@ -283,8 +291,9 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
             # need to ensure dtype=int64 in case range is empty on Python 2
             return pd.Index(range(self.sizes[key]), name=key, dtype=np.int64)
 
-    def _calc_assign_results(self, kwargs):
-        results = SortedKeysDict()
+    def _calc_assign_results(self, kwargs: Mapping[str, T]
+                             ) -> MutableMapping[str, T]:
+        results = SortedKeysDict()  # type: SortedKeysDict[str, T]
         for k, v in kwargs.items():
             if callable(v):
                 results[k] = v(self)
@@ -372,7 +381,8 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
         out.attrs.update(*args, **kwargs)
         return out
 
-    def pipe(self, func, *args, **kwargs):
+    def pipe(self, func: Union[Callable[..., T], Tuple[Callable[..., T], str]],
+             *args, **kwargs) -> T:
         """
         Apply func(self, *args, **kwargs)
 
@@ -424,15 +434,14 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
         if isinstance(func, tuple):
             func, target = func
             if target in kwargs:
-                msg = ('%s is both the pipe target and a keyword argument'
-                       % target)
-                raise ValueError(msg)
+                raise ValueError('%s is both the pipe target and a keyword '
+                                 'argument' % target)
             kwargs[target] = self
             return func(*args, **kwargs)
         else:
             return func(self, *args, **kwargs)
 
-    def groupby(self, group, squeeze=True):
+    def groupby(self, group, squeeze: bool = True):
         """Returns a GroupBy object for performing grouped operations.
 
         Parameters
@@ -478,8 +487,9 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
         """  # noqa
         return self._groupby_cls(self, group, squeeze=squeeze)
 
-    def groupby_bins(self, group, bins, right=True, labels=None, precision=3,
-                     include_lowest=False, squeeze=True):
+    def groupby_bins(self, group, bins, right: bool = True, labels=None,
+                     precision: int = 3, include_lowest: bool = False,
+                     squeeze: bool = True):
         """Returns a GroupBy object for performing grouped operations.
 
         Rather than using all unique values of `group`, the values are discretized
@@ -530,7 +540,9 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
                                              'precision': precision,
                                              'include_lowest': include_lowest})
 
-    def rolling(self, dim=None, min_periods=None, center=False, **dim_kwargs):
+    def rolling(self, dim: Optional[Mapping[Hashable, int]] = None,
+                min_periods: Optional[int] = None, center: bool = False,
+                **dim_kwargs: int):
         """
         Rolling window object.
 
@@ -590,8 +602,11 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
         return self._rolling_cls(self, dim, min_periods=min_periods,
                                  center=center)
 
-    def coarsen(self, dim=None, boundary='exact', side='left',
-                coord_func='mean', **dim_kwargs):
+    def coarsen(self, dim: Optional[Mapping[Hashable, int]] = None,
+                boundary: str = 'exact',
+                side: Union[str, Mapping[Hashable, str]] = 'left',
+                coord_func: str = 'mean',
+                **dim_kwargs: int):
         """
         Coarsen object.
 
@@ -650,8 +665,12 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
             self, dim, boundary=boundary, side=side,
             coord_func=coord_func)
 
-    def resample(self, indexer=None, skipna=None, closed=None, label=None,
-                 base=0, keep_attrs=None, loffset=None, **indexer_kwargs):
+    def resample(self, indexer: Optional[Mapping[Hashable, str]] = None,
+                 skipna=None, closed: Optional[str] = None,
+                 label: Optional[str] = None,
+                 base: int = 0, keep_attrs: Optional[bool] = None,
+                 loffset=None,
+                 **indexer_kwargs: str):
         """Returns a Resample object for performing resampling operations.
 
         Handles both downsampling and upsampling. If any intervals contain no
@@ -745,12 +764,11 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
                 "objects, e.g., data.resample(time='1D').mean()")
 
         indexer = either_dict_or_kwargs(indexer, indexer_kwargs, 'resample')
-
         if len(indexer) != 1:
             raise ValueError(
                 "Resampling only supported along single dimensions."
             )
-        dim, freq = indexer.popitem()
+        dim, freq = next(iter(indexer.items()))
 
         dim_name = dim
         dim_coord = self[dim]
@@ -772,7 +790,7 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
 
         return resampler
 
-    def where(self, cond, other=dtypes.NA, drop=False):
+    def where(self, cond, other=dtypes.NA, drop: bool = False):
         """Filter elements from this object according to a condition.
 
         This operation follows the normal broadcasting and alignment rules that
@@ -858,7 +876,7 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
 
         return ops.where_method(self, cond, other)
 
-    def close(self):
+    def close(self: Any) -> None:
         """Close any files linked to this object
         """
         if self._file_obj is not None:
@@ -914,14 +932,14 @@ class DataWithCoords(SupportsArithmetic, AttrAccessMixin):
             dask='allowed',
         )
 
-    def __enter__(self):
+    def __enter__(self: T) -> T:
         return self
 
-    def __exit__(self, exc_type, exc_value, traceback):
+    def __exit__(self, exc_type, exc_value, traceback) -> None:
         self.close()
 
 
-def full_like(other, fill_value, dtype=None):
+def full_like(other, fill_value, dtype: Union[str, np.dtype, None] = None):
     """Return a new object with the same shape and type as a given object.
 
     Parameters
@@ -961,7 +979,8 @@ def full_like(other, fill_value, dtype=None):
         raise TypeError("Expected DataArray, Dataset, or Variable")
 
 
-def _full_like_variable(other, fill_value, dtype=None):
+def _full_like_variable(other, fill_value,
+                        dtype: Union[str, np.dtype, None] = None):
     """Inner function of full_like, where other must be a variable
     """
     from .variable import Variable
@@ -978,27 +997,28 @@ def _full_like_variable(other, fill_value, dtype=None):
     return Variable(dims=other.dims, data=data, attrs=other.attrs)
 
 
-def zeros_like(other, dtype=None):
+def zeros_like(other, dtype: Union[str, np.dtype, None] = None):
     """Shorthand for full_like(other, 0, dtype)
     """
     return full_like(other, 0, dtype)
 
 
-def ones_like(other, dtype=None):
+def ones_like(other, dtype: Union[str, np.dtype, None] = None):
     """Shorthand for full_like(other, 1, dtype)
     """
     return full_like(other, 1, dtype)
 
 
-def is_np_datetime_like(dtype):
+def is_np_datetime_like(dtype: Union[str, np.dtype]) -> bool:
     """Check if a dtype is a subclass of the numpy datetime types
     """
     return (np.issubdtype(dtype, np.datetime64) or
             np.issubdtype(dtype, np.timedelta64))
 
 
-def _contains_cftime_datetimes(array):
-    """Check if an array contains cftime.datetime objects"""
+def _contains_cftime_datetimes(array) -> bool:
+    """Check if an array contains cftime.datetime objects
+    """
     try:
         from cftime import datetime as cftime_datetime
     except ImportError:
@@ -1015,12 +1035,14 @@ def _contains_cftime_datetimes(array):
             return False
 
 
-def contains_cftime_datetimes(var):
-    """Check if an xarray.Variable contains cftime.datetime objects"""
+def contains_cftime_datetimes(var) -> bool:
+    """Check if an xarray.Variable contains cftime.datetime objects
+    """
     return _contains_cftime_datetimes(var.data)
 
 
-def _contains_datetime_like_objects(var):
+def _contains_datetime_like_objects(var) -> bool:
     """Check if a variable contains datetime like objects (either
-    np.datetime64, np.timedelta64, or cftime.datetime)"""
+    np.datetime64, np.timedelta64, or cftime.datetime)
+    """
     return is_np_datetime_like(var.dtype) or contains_cftime_datetimes(var)
diff --git a/xarray/core/dataarray.py b/xarray/core/dataarray.py
index c24703f5..82e20d12 100644
--- a/xarray/core/dataarray.py
+++ b/xarray/core/dataarray.py
@@ -231,9 +231,6 @@ class DataArray(AbstractArray, DataWithCoords):
             coords, dims = _infer_coords_and_dims(data.shape, coords, dims)
             variable = Variable(dims, data, attrs, encoding, fastpath=True)
 
-        # uncomment for a useful consistency check:
-        # assert all(isinstance(v, Variable) for v in coords.values())
-
         # These fully describe a DataArray
         self._variable = variable
         self._coords = coords
@@ -1402,7 +1399,7 @@ class DataArray(AbstractArray, DataWithCoords):
         ds = self._to_temp_dataset().unstack(dim)
         return self._from_temp_dataset(ds)
 
-    def transpose(self, *dims):
+    def transpose(self, *dims) -> 'DataArray':
         """Return a new DataArray object with transposed dimensions.
 
         Parameters
@@ -1430,6 +1427,10 @@ class DataArray(AbstractArray, DataWithCoords):
         variable = self.variable.transpose(*dims)
         return self._replace(variable)
 
+    @property
+    def T(self) -> 'DataArray':
+        return self.transpose()
+
     def drop(self, labels, dim=None):
         """Drop coordinates or index labels from this DataArray.
 
diff --git a/xarray/core/dataset.py b/xarray/core/dataset.py
index e3b2e3c3..cf6631fa 100644
--- a/xarray/core/dataset.py
+++ b/xarray/core/dataset.py
@@ -938,6 +938,7 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
         """
         variables = OrderedDict()  # type: OrderedDict[Any, Variable]
         coord_names = set()
+        indexes = OrderedDict()  # type: OrderedDict[Any, pd.Index]
 
         for name in names:
             try:
@@ -948,6 +949,8 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
                 variables[var_name] = var
                 if ref_name in self._coord_names or ref_name in self.dims:
                     coord_names.add(var_name)
+                if (var_name,) == var.dims:
+                    indexes[var_name] = var.to_index()
 
         needed_dims = set()  # type: set
         for v in variables.values():
@@ -959,12 +962,8 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
             if set(self.variables[k].dims) <= needed_dims:
                 variables[k] = self._variables[k]
                 coord_names.add(k)
-
-        if self._indexes is None:
-            indexes = None
-        else:
-            indexes = OrderedDict((k, v) for k, v in self._indexes.items()
-                                  if k in coord_names)
+                if k in self.indexes:
+                    indexes[k] = self.indexes[k]
 
         return self._replace(variables, coord_names, dims, indexes=indexes)
 
@@ -1503,9 +1502,13 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
             raise ValueError("dimensions %r do not exist" % invalid)
 
         # all indexers should be int, slice, np.ndarrays, or Variable
-        indexers_list = []
+        indexers_list = []  # type: List[Tuple[Any, Union[slice, Variable]]]
         for k, v in indexers.items():
-            if isinstance(v, (slice, Variable)):
+            if isinstance(v, slice):
+                indexers_list.append((k, v))
+                continue
+
+            if isinstance(v, Variable):
                 pass
             elif isinstance(v, DataArray):
                 v = v.variable
@@ -1524,14 +1527,19 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
                         v = _parse_array_of_cftime_strings(v, index.date_type)
 
                 if v.ndim == 0:
-                    v = as_variable(v)
+                    v = Variable((), v)
                 elif v.ndim == 1:
-                    v = as_variable((k, v))
+                    v = IndexVariable((k,), v)
                 else:
                     raise IndexError(
                         "Unlabeled multi-dimensional array cannot be "
                         "used for indexing: {}".format(k))
+
+            if v.ndim == 1:
+                v = v.to_index_variable()
+
             indexers_list.append((k, v))
+
         return indexers_list
 
     def _get_indexers_coords_and_indexes(self, indexers):
@@ -1631,7 +1639,7 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
 
             if name in self.indexes:
                 new_var, new_index = isel_variable_and_index(
-                    var, self.indexes[name], var_indexers)
+                    name, var, self.indexes[name], var_indexers)
                 if new_index is not None:
                     indexes[name] = new_index
             else:
@@ -2117,15 +2125,20 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
         indexes = OrderedDict(
             (k, v) for k, v in obj.indexes.items() if k not in indexers)
         selected = self._replace_with_new_dims(
-            variables, coord_names, indexes=indexes)
+            variables.copy(), coord_names, indexes=indexes)
 
         # attach indexer as coordinate
         variables.update(indexers)
+        indexes.update(
+            (k, v.to_index()) for k, v in indexers.items() if v.dims == (k,)
+        )
+
         # Extract coordinates from indexers
         coord_vars, new_indexes = (
             selected._get_indexers_coords_and_indexes(coords))
         variables.update(coord_vars)
         indexes.update(new_indexes)
+
         coord_names = (set(variables)
                        .intersection(obj._coord_names)
                        .union(coord_vars))
@@ -2401,6 +2414,7 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
                     ' variable name.'.format(dim=d))
 
         variables = OrderedDict()
+        coord_names = self._coord_names.copy()
         # If dim is a dict, then ensure that the values are either integers
         # or iterables.
         for k, v in dim.items():
@@ -2410,7 +2424,7 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
                 # value within the dim dict to the length of the iterable
                 # for later use.
                 variables[k] = xr.IndexVariable((k,), v)
-                self._coord_names.add(k)
+                coord_names.add(k)
                 dim[k] = variables[k].size
             elif isinstance(v, int):
                 pass  # Do nothing if the dimensions value is just an int
@@ -2420,7 +2434,7 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
 
         for k, v in self._variables.items():
             if k not in dim:
-                if k in self._coord_names:  # Do not change coordinates
+                if k in coord_names:  # Do not change coordinates
                     variables[k] = v
                 else:
                     result_ndim = len(v.dims) + len(axis)
@@ -2452,10 +2466,10 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
                 variables[k] = v.set_dims(k)
 
         new_dims = self._dims.copy()
-        for d in dim:
-            new_dims[d] = 1
+        new_dims.update(dim)
 
-        return self._replace(variables, dims=new_dims)
+        return self._replace_vars_and_dims(
+            variables, dims=new_dims, coord_names=coord_names)
 
     def set_index(self, indexes=None, append=False, inplace=None,
                   **indexes_kwargs):
diff --git a/xarray/core/indexes.py b/xarray/core/indexes.py
index 6d8b5530..eccb72b6 100644
--- a/xarray/core/indexes.py
+++ b/xarray/core/indexes.py
@@ -1,6 +1,6 @@
 import collections.abc
 from collections import OrderedDict
-from typing import Any, Iterable, Mapping, Optional, Tuple, Union
+from typing import Any, Hashable, Iterable, Mapping, Optional, Tuple, Union
 
 import pandas as pd
 
@@ -59,6 +59,7 @@ def default_indexes(
 
 
 def isel_variable_and_index(
+    name: Hashable,
     variable: Variable,
     index: pd.Index,
     indexers: Mapping[Any, Union[slice, Variable]],
@@ -75,8 +76,8 @@ def isel_variable_and_index(
 
     new_variable = variable.isel(indexers)
 
-    if new_variable.ndim != 1:
-        # can't preserve a index if result is not 0D
+    if new_variable.dims != (name,):
+        # can't preserve a index if result has new dimensions
         return new_variable, None
 
     # we need to compute the new index
diff --git a/xarray/core/pycompat.py b/xarray/core/pycompat.py
index 0df0e727..aa2cc5a0 100644
--- a/xarray/core/pycompat.py
+++ b/xarray/core/pycompat.py
@@ -1,6 +1,5 @@
-# flake8: noqa
 import sys
-import typing
+from collections import abc
 
 import numpy as np
 
@@ -13,6 +12,38 @@ try:
 except ImportError:  # pragma: no cover
     dask_array_type = ()
 
-# Ensure we have some more recent additions to the typing module.
-# Note that TYPE_CHECKING itself is not available on Python 3.5.1.
-TYPE_CHECKING = sys.version >= '3.5.3' and typing.TYPE_CHECKING
+
+if sys.version < '3.5.3':
+    TYPE_CHECKING = False
+
+    class _ABCDummyBrackets(type(abc.Mapping)):  # abc.ABCMeta
+        def __getitem__(cls, name):
+            return cls
+
+    class Mapping(abc.Mapping, metaclass=_ABCDummyBrackets):
+        pass
+
+    class MutableMapping(abc.MutableMapping, metaclass=_ABCDummyBrackets):
+        pass
+
+    class MutableSet(abc.MutableSet, metaclass=_ABCDummyBrackets):
+        pass
+
+else:
+    from typing import TYPE_CHECKING  # noqa: F401
+
+    # from typing import Mapping, MutableMapping, MutableSet
+
+    # The above confuses mypy 0.700;
+    # see: https://github.com/python/mypy/issues/6652
+    # As a workaround, use:
+    #
+    # from typing import Mapping, MutableMapping, MutableSet
+    # try:
+    #     from .pycompat import Mapping, MutableMapping, MutableSet
+    # except ImportError:
+    #      pass
+    #
+    # This is only necessary in modules that define subclasses of the
+    # abstract collections; when only type inference is needed, one can just
+    # use typing also in Python 3.5.0~3.5.2 (although mypy will misbehave).
diff --git a/xarray/core/utils.py b/xarray/core/utils.py
index 349c8f98..94787dd3 100644
--- a/xarray/core/utils.py
+++ b/xarray/core/utils.py
@@ -7,15 +7,28 @@ import os.path
 import re
 import warnings
 from collections import OrderedDict
-from collections.abc import Iterable, Mapping, MutableMapping, MutableSet
+from typing import (AbstractSet, Any, Callable, Container, Dict, Hashable,
+                    Iterable, Iterator, Optional, Sequence,
+                    Tuple, TypeVar, cast)
 
 import numpy as np
 import pandas as pd
 
 from .pycompat import dask_array_type
 
+from typing import Mapping, MutableMapping, MutableSet
+try:  # Fix typed collections in Python 3.5.0~3.5.2
+    from .pycompat import Mapping, MutableMapping, MutableSet  # noqa: F811
+except ImportError:
+    pass
 
-def _check_inplace(inplace, default=False):
+
+K = TypeVar('K')
+V = TypeVar('V')
+T = TypeVar('T')
+
+
+def _check_inplace(inplace: bool, default: bool = False) -> bool:
     if inplace is None:
         inplace = default
     else:
@@ -26,16 +39,16 @@ def _check_inplace(inplace, default=False):
     return inplace
 
 
-def alias_message(old_name, new_name):
+def alias_message(old_name: str, new_name: str) -> str:
     return '%s has been deprecated. Use %s instead.' % (old_name, new_name)
 
 
-def alias_warning(old_name, new_name, stacklevel=3):
+def alias_warning(old_name: str, new_name: str, stacklevel: int = 3) -> None:
     warnings.warn(alias_message(old_name, new_name), FutureWarning,
                   stacklevel=stacklevel)
 
 
-def alias(obj, old_name):
+def alias(obj: Callable[..., T], old_name: str) -> Callable[..., T]:
     assert isinstance(old_name, str)
 
     @functools.wraps(obj)
@@ -46,7 +59,7 @@ def alias(obj, old_name):
     return wrapper
 
 
-def _maybe_cast_to_cftimeindex(index):
+def _maybe_cast_to_cftimeindex(index: pd.Index) -> pd.Index:
     from ..coding.cftimeindex import CFTimeIndex
 
     if index.dtype == 'O':
@@ -58,7 +71,7 @@ def _maybe_cast_to_cftimeindex(index):
         return index
 
 
-def safe_cast_to_index(array):
+def safe_cast_to_index(array: Any) -> pd.Index:
     """Given an array, safely cast it to a pandas.Index.
 
     If it is already a pandas.Index, return it unchanged.
@@ -79,7 +92,9 @@ def safe_cast_to_index(array):
     return _maybe_cast_to_cftimeindex(index)
 
 
-def multiindex_from_product_levels(levels, names=None):
+def multiindex_from_product_levels(levels: Sequence[pd.Index],
+                                   names: Optional[Sequence[str]] = None
+                                   ) -> pd.MultiIndex:
     """Creating a MultiIndex from a product without refactorizing levels.
 
     Keeping levels the same gives back the original labels when we unstack.
@@ -117,7 +132,7 @@ def maybe_wrap_array(original, new_array):
         return new_array
 
 
-def equivalent(first, second):
+def equivalent(first: T, second: T) -> bool:
     """Compare two objects for equivalence (identity or equality), using
     array_equiv if either object is an ndarray
     """
@@ -131,8 +146,8 @@ def equivalent(first, second):
                 (pd.isnull(first) and pd.isnull(second)))
 
 
-def peek_at(iterable):
-    """Returns the first value from iterable, as well as a new iterable with
+def peek_at(iterable: Iterable[T]) -> Tuple[T, Iterator[T]]:
+    """Returns the first value from iterable, as well as a new iterator with
     the same content as the original iterable
     """
     gen = iter(iterable)
@@ -140,7 +155,9 @@ def peek_at(iterable):
     return peek, itertools.chain([peek], gen)
 
 
-def update_safety_check(first_dict, second_dict, compat=equivalent):
+def update_safety_check(first_dict: MutableMapping[K, V],
+                        second_dict: Mapping[K, V],
+                        compat: Callable[[V, V], bool] = equivalent) -> None:
     """Check the safety of updating one dictionary with another.
 
     Raises ValueError if dictionaries have non-compatible values for any key,
@@ -162,7 +179,10 @@ def update_safety_check(first_dict, second_dict, compat=equivalent):
                              'overriding values; conflicting key %r' % k)
 
 
-def remove_incompatible_items(first_dict, second_dict, compat=equivalent):
+def remove_incompatible_items(first_dict: MutableMapping[K, V],
+                              second_dict: Mapping[K, V],
+                              compat: Callable[[V, V], bool] = equivalent
+                              ) -> None:
     """Remove incompatible items from the first dictionary in-place.
 
     Items are retained if their keys are found in both dictionaries and the
@@ -177,21 +197,22 @@ def remove_incompatible_items(first_dict, second_dict, compat=equivalent):
         checks for equivalence.
     """
     for k in list(first_dict):
-        if (k not in second_dict or
-            (k in second_dict and
-                not compat(first_dict[k], second_dict[k]))):
+        if k not in second_dict or not compat(first_dict[k], second_dict[k]):
             del first_dict[k]
 
 
-def is_dict_like(value):
+def is_dict_like(value: Any) -> bool:
     return hasattr(value, 'keys') and hasattr(value, '__getitem__')
 
 
-def is_full_slice(value):
+def is_full_slice(value: Any) -> bool:
     return isinstance(value, slice) and value == slice(None)
 
 
-def either_dict_or_kwargs(pos_kwargs, kw_kwargs, func_name):
+def either_dict_or_kwargs(pos_kwargs: Optional[Mapping[Hashable, T]],
+                          kw_kwargs: Mapping[str, T],
+                          func_name: str
+                          ) -> Mapping[Hashable, T]:
     if pos_kwargs is not None:
         if not is_dict_like(pos_kwargs):
             raise ValueError('the first argument to .%s must be a dictionary'
@@ -201,10 +222,12 @@ def either_dict_or_kwargs(pos_kwargs, kw_kwargs, func_name):
                              'arguments to .%s' % func_name)
         return pos_kwargs
     else:
-        return kw_kwargs
+        # Need an explicit cast to appease mypy due to invariance; see
+        # https://github.com/python/mypy/issues/6228
+        return cast(Mapping[Hashable, T], kw_kwargs)
 
 
-def is_scalar(value):
+def is_scalar(value: Any) -> bool:
     """Whether to treat a value as a scalar.
 
     Any non-iterable, string, or 0-D array
@@ -215,7 +238,7 @@ def is_scalar(value):
         isinstance(value, (Iterable, ) + dask_array_type))
 
 
-def is_valid_numpy_dtype(dtype):
+def is_valid_numpy_dtype(dtype: Any) -> bool:
     try:
         np.dtype(dtype)
     except (TypeError, ValueError):
@@ -224,15 +247,17 @@ def is_valid_numpy_dtype(dtype):
         return True
 
 
-def to_0d_object_array(value):
-    """Given a value, wrap it in a 0-D numpy.ndarray with dtype=object."""
+def to_0d_object_array(value: Any) -> np.ndarray:
+    """Given a value, wrap it in a 0-D numpy.ndarray with dtype=object.
+    """
     result = np.empty((), dtype=object)
     result[()] = value
     return result
 
 
-def to_0d_array(value):
-    """Given a value, wrap it in a 0-D numpy.ndarray."""
+def to_0d_array(value: Any) -> np.ndarray:
+    """Given a value, wrap it in a 0-D numpy.ndarray.
+    """
     if np.isscalar(value) or (isinstance(value, np.ndarray) and
                               value.ndim == 0):
         return np.array(value)
@@ -240,7 +265,8 @@ def to_0d_array(value):
         return to_0d_object_array(value)
 
 
-def dict_equiv(first, second, compat=equivalent):
+def dict_equiv(first: Mapping[K, V], second: Mapping[K, V],
+               compat: Callable[[V, V], bool] = equivalent) -> bool:
     """Test equivalence of two dict-like objects. If any of the values are
     numpy arrays, compare them correctly.
 
@@ -266,7 +292,10 @@ def dict_equiv(first, second, compat=equivalent):
     return True
 
 
-def ordered_dict_intersection(first_dict, second_dict, compat=equivalent):
+def ordered_dict_intersection(first_dict: Mapping[K, V],
+                              second_dict: Mapping[K, V],
+                              compat: Callable[[V, V], bool] = equivalent
+                              ) -> MutableMapping[K, V]:
     """Return the intersection of two dictionaries as a new OrderedDict.
 
     Items are retained if their keys are found in both dictionaries and the
@@ -290,11 +319,10 @@ def ordered_dict_intersection(first_dict, second_dict, compat=equivalent):
     return new_dict
 
 
-class SingleSlotPickleMixin(object):
+class SingleSlotPickleMixin:
     """Mixin class to add the ability to pickle objects whose state is defined
     by a single __slots__ attribute. Only necessary under Python 2.
     """
-
     def __getstate__(self):
         return getattr(self, self.__slots__[0])
 
@@ -302,160 +330,123 @@ class SingleSlotPickleMixin(object):
         setattr(self, self.__slots__[0], state)
 
 
-class Frozen(Mapping, SingleSlotPickleMixin):
+class Frozen(Mapping[K, V], SingleSlotPickleMixin):
     """Wrapper around an object implementing the mapping interface to make it
     immutable. If you really want to modify the mapping, the mutable version is
     saved under the `mapping` attribute.
     """
     __slots__ = ['mapping']
 
-    def __init__(self, mapping):
+    def __init__(self, mapping: Mapping[K, V]):
         self.mapping = mapping
 
-    def __getitem__(self, key):
+    def __getitem__(self, key: K) -> V:
         return self.mapping[key]
 
-    def __iter__(self):
+    def __iter__(self) -> Iterator[K]:
         return iter(self.mapping)
 
-    def __len__(self):
+    def __len__(self) -> int:
         return len(self.mapping)
 
-    def __contains__(self, key):
+    def __contains__(self, key: object) -> bool:
         return key in self.mapping
 
-    def __repr__(self):
+    def __repr__(self) -> str:
         return '%s(%r)' % (type(self).__name__, self.mapping)
 
 
-def FrozenOrderedDict(*args, **kwargs):
+def FrozenOrderedDict(*args, **kwargs) -> Frozen:
     return Frozen(OrderedDict(*args, **kwargs))
 
 
-class SortedKeysDict(MutableMapping, SingleSlotPickleMixin):
+class SortedKeysDict(MutableMapping[K, V], SingleSlotPickleMixin):
     """An wrapper for dictionary-like objects that always iterates over its
     items in sorted order by key but is otherwise equivalent to the underlying
     mapping.
     """
     __slots__ = ['mapping']
 
-    def __init__(self, mapping=None):
+    def __init__(self, mapping: Optional[MutableMapping[K, V]] = None):
         self.mapping = {} if mapping is None else mapping
 
-    def __getitem__(self, key):
+    def __getitem__(self, key: K) -> V:
         return self.mapping[key]
 
-    def __setitem__(self, key, value):
+    def __setitem__(self, key: K, value: V) -> None:
         self.mapping[key] = value
 
-    def __delitem__(self, key):
+    def __delitem__(self, key: K) -> None:
         del self.mapping[key]
 
-    def __iter__(self):
+    def __iter__(self) -> Iterator[K]:
         return iter(sorted(self.mapping))
 
-    def __len__(self):
+    def __len__(self) -> int:
         return len(self.mapping)
 
-    def __contains__(self, key):
+    def __contains__(self, key: object) -> bool:
         return key in self.mapping
 
-    def __repr__(self):
+    def __repr__(self) -> str:
         return '%s(%r)' % (type(self).__name__, self.mapping)
 
-    def copy(self):
-        return type(self)(self.mapping.copy())
-
-
-class ChainMap(MutableMapping, SingleSlotPickleMixin):
-    """Partial backport of collections.ChainMap from Python>=3.3
-
-    Don't return this from any public APIs, since some of the public methods
-    for a MutableMapping are missing (they will raise a NotImplementedError)
-    """
-    __slots__ = ['maps']
-
-    def __init__(self, *maps):
-        self.maps = maps
-
-    def __getitem__(self, key):
-        for mapping in self.maps:
-            try:
-                return mapping[key]
-            except KeyError:
-                pass
-        raise KeyError(key)
 
-    def __setitem__(self, key, value):
-        self.maps[0][key] = value
-
-    def __delitem__(self, value):  # pragma: no cover
-        raise NotImplementedError
-
-    def __iter__(self):
-        seen = set()
-        for mapping in self.maps:
-            for item in mapping:
-                if item not in seen:
-                    yield item
-                    seen.add(item)
-
-    def __len__(self):
-        raise len(iter(self))
-
-
-class OrderedSet(MutableSet):
+class OrderedSet(MutableSet[T]):
     """A simple ordered set.
 
     The API matches the builtin set, but it preserves insertion order of
     elements, like an OrderedDict.
     """
-
-    def __init__(self, values=None):
-        self._ordered_dict = OrderedDict()
+    def __init__(self, values: Optional[AbstractSet[T]] = None):
+        self._ordered_dict = OrderedDict()  # type: MutableMapping[T, None]
         if values is not None:
-            self |= values
+            # Disable type checking - both mypy and PyCharm believes that
+            # we're altering the type of self in place (see signature of
+            # MutableSet.__ior__)
+            self |= values  # type: ignore
 
     # Required methods for MutableSet
 
-    def __contains__(self, value):
+    def __contains__(self, value: object) -> bool:
         return value in self._ordered_dict
 
-    def __iter__(self):
+    def __iter__(self) -> Iterator[T]:
         return iter(self._ordered_dict)
 
-    def __len__(self):
+    def __len__(self) -> int:
         return len(self._ordered_dict)
 
-    def add(self, value):
+    def add(self, value: T) -> None:
         self._ordered_dict[value] = None
 
-    def discard(self, value):
+    def discard(self, value: T) -> None:
         del self._ordered_dict[value]
 
     # Additional methods
 
-    def update(self, values):
-        self |= values
+    def update(self, values: AbstractSet[T]) -> None:
+        # See comment on __init__ re. type checking
+        self |= values  # type: ignore
 
-    def __repr__(self):
+    def __repr__(self) -> str:
         return '%s(%r)' % (type(self).__name__, list(self))
 
 
-class NdimSizeLenMixin(object):
+class NdimSizeLenMixin:
     """Mixin class that extends a class that defines a ``shape`` property to
     one that also defines ``ndim``, ``size`` and ``__len__``.
     """
     @property
-    def ndim(self):
+    def ndim(self: Any) -> int:
         return len(self.shape)
 
     @property
-    def size(self):
+    def size(self: Any) -> int:
         # cast to int so that shape = () gives size = 1
         return int(np.prod(self.shape))
 
-    def __len__(self):
+    def __len__(self: Any) -> int:
         try:
             return self.shape[0]
         except IndexError:
@@ -470,27 +461,27 @@ class NDArrayMixin(NdimSizeLenMixin):
     `dtype`, `shape` and `__getitem__`.
     """
     @property
-    def dtype(self):
+    def dtype(self: Any) -> np.dtype:
         return self.array.dtype
 
     @property
-    def shape(self):
+    def shape(self: Any) -> Tuple[int]:
         return self.array.shape
 
-    def __getitem__(self, key):
+    def __getitem__(self: Any, key):
         return self.array[key]
 
-    def __repr__(self):
+    def __repr__(self: Any) -> str:
         return '%s(array=%r)' % (type(self).__name__, self.array)
 
 
-class ReprObject(object):
-    """Object that prints as the given value, for use with sentinel values."""
-
+class ReprObject:
+    """Object that prints as the given value, for use with sentinel values.
+    """
     def __init__(self, value: str):
         self._value = value
 
-    def __repr__(self):
+    def __repr__(self) -> str:
         return self._value
 
 
@@ -506,16 +497,16 @@ def close_on_error(f):
         raise
 
 
-def is_remote_uri(path):
+def is_remote_uri(path: str) -> bool:
     return bool(re.search(r'^https?\://', path))
 
 
-def is_grib_path(path):
+def is_grib_path(path: str) -> bool:
     _, ext = os.path.splitext(path)
     return ext in ['.grib', '.grb', '.grib2', '.grb2']
 
 
-def is_uniform_spaced(arr, **kwargs):
+def is_uniform_spaced(arr, **kwargs) -> bool:
     """Return True if values of an array are uniformly spaced and sorted.
 
     >>> is_uniform_spaced(range(5))
@@ -527,11 +518,12 @@ def is_uniform_spaced(arr, **kwargs):
     """
     arr = np.array(arr, dtype=float)
     diffs = np.diff(arr)
-    return np.isclose(diffs.min(), diffs.max(), **kwargs)
+    return bool(np.isclose(diffs.min(), diffs.max(), **kwargs))
 
 
-def hashable(v):
-    """Determine whether `v` can be hashed."""
+def hashable(v: Any) -> bool:
+    """Determine whether `v` can be hashed.
+    """
     try:
         hash(v)
     except TypeError:
@@ -543,9 +535,10 @@ def not_implemented(*args, **kwargs):
     return NotImplemented
 
 
-def decode_numpy_dict_values(attrs):
+def decode_numpy_dict_values(attrs: Mapping[K, V]) -> Dict[K, V]:
     """Convert attribute values from numpy objects to native Python objects,
-    for use in to_dict"""
+    for use in to_dict
+    """
     attrs = dict(attrs)
     for k, v in attrs.items():
         if isinstance(v, np.ndarray):
@@ -565,46 +558,43 @@ def ensure_us_time_resolution(val):
     return val
 
 
-class HiddenKeyDict(MutableMapping):
-    '''
-    Acts like a normal dictionary, but hides certain keys.
-    '''
+class HiddenKeyDict(MutableMapping[K, V]):
+    """Acts like a normal dictionary, but hides certain keys.
+    """
     # ``__init__`` method required to create instance from class.
 
-    def __init__(self, data, hidden_keys):
+    def __init__(self, data: MutableMapping[K, V], hidden_keys: Iterable[K]):
         self._data = data
-        if type(hidden_keys) not in (list, tuple):
-            raise TypeError("hidden_keys must be a list or tuple")
-        self._hidden_keys = hidden_keys
+        self._hidden_keys = frozenset(hidden_keys)
 
-    def _raise_if_hidden(self, key):
+    def _raise_if_hidden(self, key: K) -> None:
         if key in self._hidden_keys:
             raise KeyError('Key `%r` is hidden.' % key)
 
     # The next five methods are requirements of the ABC.
-    def __setitem__(self, key, value):
+    def __setitem__(self, key: K, value: V) -> None:
         self._raise_if_hidden(key)
         self._data[key] = value
 
-    def __getitem__(self, key):
+    def __getitem__(self, key: K) -> V:
         self._raise_if_hidden(key)
         return self._data[key]
 
-    def __delitem__(self, key):
+    def __delitem__(self, key: K) -> None:
         self._raise_if_hidden(key)
         del self._data[key]
 
-    def __iter__(self):
+    def __iter__(self) -> Iterator[K]:
         for k in self._data:
             if k not in self._hidden_keys:
                 yield k
 
-    def __len__(self):
-        num_hidden = sum([k in self._hidden_keys for k in self._data])
+    def __len__(self) -> int:
+        num_hidden = len(self._hidden_keys & self._data.keys())
         return len(self._data) - num_hidden
 
 
-def get_temp_dimname(dims, new_dim):
+def get_temp_dimname(dims: Container[Hashable], new_dim: Hashable) -> Hashable:
     """ Get an new dimension name based on new_dim, that is not used in dims.
     If the same name exists, we add an underscore(s) in the head.
 
@@ -618,5 +608,5 @@ def get_temp_dimname(dims, new_dim):
         -> ['__rolling']
     """
     while new_dim in dims:
-        new_dim = '_' + new_dim
+        new_dim = '_' + str(new_dim)
     return new_dim
diff --git a/xarray/core/variable.py b/xarray/core/variable.py
index d6b64e7d..96c6b7bd 100644
--- a/xarray/core/variable.py
+++ b/xarray/core/variable.py
@@ -1121,7 +1121,7 @@ class Variable(common.AbstractArray, arithmetic.SupportsArithmetic,
             result = result._roll_one_dim(dim, count)
         return result
 
-    def transpose(self, *dims):
+    def transpose(self, *dims) -> 'Variable':
         """Return a new Variable object with transposed dimensions.
 
         Parameters
@@ -1155,6 +1155,10 @@ class Variable(common.AbstractArray, arithmetic.SupportsArithmetic,
         return type(self)(dims, data, self._attrs, self._encoding,
                           fastpath=True)
 
+    @property
+    def T(self) -> 'Variable':
+        return self.transpose()
+
     def expand_dims(self, *args):
         import warnings
         warnings.warn('Variable.expand_dims is deprecated: use '


Test patch:
diff --git a/xarray/testing.py b/xarray/testing.py
index 794c0614..eb8a0e86 100644
--- a/xarray/testing.py
+++ b/xarray/testing.py
@@ -1,8 +1,12 @@
 """Testing functions exposed to the user API"""
+from collections import OrderedDict
+
 import numpy as np
+import pandas as pd
 
 from xarray.core import duck_array_ops
 from xarray.core import formatting
+from xarray.core.indexes import default_indexes
 
 
 def _decode_string_data(data):
@@ -143,8 +147,37 @@ def assert_allclose(a, b, rtol=1e-05, atol=1e-08, decode_bytes=True):
                         .format(type(a)))
 
 
-def assert_combined_tile_ids_equal(dict1, dict2):
-    assert len(dict1) == len(dict2)
-    for k, v in dict1.items():
-        assert k in dict2.keys()
-        assert_equal(dict1[k], dict2[k])
+def _assert_indexes_invariants_checks(indexes, possible_coord_variables, dims):
+    import xarray as xr
+
+    assert isinstance(indexes, OrderedDict), indexes
+    assert all(isinstance(v, pd.Index) for v in indexes.values()), \
+        {k: type(v) for k, v in indexes.items()}
+
+    index_vars = {k for k, v in possible_coord_variables.items()
+                  if isinstance(v, xr.IndexVariable)}
+    assert indexes.keys() <= index_vars, (set(indexes), index_vars)
+
+    # Note: when we support non-default indexes, these checks should be opt-in
+    # only!
+    defaults = default_indexes(possible_coord_variables, dims)
+    assert indexes.keys() == defaults.keys(), \
+        (set(indexes), set(defaults))
+    assert all(v.equals(defaults[k]) for k, v in indexes.items()), \
+        (indexes, defaults)
+
+
+def _assert_indexes_invariants(a):
+    """Separate helper function for checking indexes invariants only."""
+    import xarray as xr
+
+    if isinstance(a, xr.DataArray):
+        if a._indexes is not None:
+            _assert_indexes_invariants_checks(a._indexes, a._coords, a.dims)
+    elif isinstance(a, xr.Dataset):
+        if a._indexes is not None:
+            _assert_indexes_invariants_checks(
+                a._indexes, a._variables, a._dims)
+    elif isinstance(a, xr.Variable):
+        # no indexes
+        pass
diff --git a/xarray/tests/__init__.py b/xarray/tests/__init__.py
index 4ebcc29a..52536070 100644
--- a/xarray/tests/__init__.py
+++ b/xarray/tests/__init__.py
@@ -13,8 +13,7 @@ import pytest
 from xarray.core import utils
 from xarray.core.options import set_options
 from xarray.core.indexing import ExplicitlyIndexed
-from xarray.testing import (assert_equal, assert_identical,  # noqa: F401
-                            assert_allclose, assert_combined_tile_ids_equal)
+import xarray.testing
 from xarray.plot.utils import import_seaborn
 
 try:
@@ -180,3 +179,25 @@ def source_ndarray(array):
     if base is None:
         base = array
     return base
+
+
+# Internal versions of xarray's test functions that validate additional
+# invariants
+# TODO: add more invariant checks.
+
+def assert_equal(a, b):
+    xarray.testing.assert_equal(a, b)
+    xarray.testing._assert_indexes_invariants(a)
+    xarray.testing._assert_indexes_invariants(b)
+
+
+def assert_identical(a, b):
+    xarray.testing.assert_identical(a, b)
+    xarray.testing._assert_indexes_invariants(a)
+    xarray.testing._assert_indexes_invariants(b)
+
+
+def assert_allclose(a, b, **kwargs):
+    xarray.testing.assert_allclose(a, b, **kwargs)
+    xarray.testing._assert_indexes_invariants(a)
+    xarray.testing._assert_indexes_invariants(b)
diff --git a/xarray/tests/test_backends.py b/xarray/tests/test_backends.py
index a20ba2df..f76d1927 100644
--- a/xarray/tests/test_backends.py
+++ b/xarray/tests/test_backends.py
@@ -3350,6 +3350,33 @@ class TestRasterio(object):
                         assert actual_shape == expected_shape
                         assert expected_val.all() == actual_val.all()
 
+    def test_rasterio_vrt_with_transform_and_size(self):
+        # Test open_rasterio() support of WarpedVRT with transform, width and
+        # height (issue #2864)
+        import rasterio
+        from rasterio.warp import calculate_default_transform
+        from affine import Affine
+        with create_tmp_geotiff() as (tmp_file, expected):
+            with rasterio.open(tmp_file) as src:
+                # Estimate the transform, width and height
+                # for a change of resolution
+                # tmp_file initial res is (1000,2000) (default values)
+                trans, w, h = calculate_default_transform(
+                    src.crs, src.crs, src.width, src.height,
+                    resolution=500, *src.bounds)
+                with rasterio.vrt.WarpedVRT(src, transform=trans,
+                                            width=w, height=h) as vrt:
+                    expected_shape = (vrt.width, vrt.height)
+                    expected_res = vrt.res
+                    expected_transform = vrt.transform
+                    with xr.open_rasterio(vrt) as da:
+                        actual_shape = (da.sizes['x'], da.sizes['y'])
+                        actual_res = da.res
+                        actual_transform = Affine(*da.transform)
+                        assert actual_res == expected_res
+                        assert actual_shape == expected_shape
+                        assert actual_transform == expected_transform
+
     @network
     def test_rasterio_vrt_network(self):
         import rasterio
diff --git a/xarray/tests/test_cftimeindex_resample.py b/xarray/tests/test_cftimeindex_resample.py
index 636f9ef7..7aca4492 100644
--- a/xarray/tests/test_cftimeindex_resample.py
+++ b/xarray/tests/test_cftimeindex_resample.py
@@ -6,31 +6,47 @@ import pandas as pd
 import xarray as xr
 from xarray.core.resample_cftime import CFTimeGrouper
 
+
 pytest.importorskip('cftime')
 pytest.importorskip('pandas', minversion='0.24')
 
 
-@pytest.fixture(
-    params=[
-        dict(start='2004-01-01T12:07:01', periods=91, freq='3D'),
-        dict(start='1892-01-03T12:07:01', periods=15, freq='41987T'),
-        dict(start='2004-01-01T12:07:01', periods=7, freq='3Q-AUG'),
-        dict(start='1892-01-03T12:07:01', periods=10, freq='3AS-JUN')
-    ],
-    ids=['3D', '41987T', '3Q_AUG', '3AS_JUN']
-)
-def time_range_kwargs(request):
-    return request.param
-
-
-@pytest.fixture()
-def datetime_index(time_range_kwargs):
-    return pd.date_range(**time_range_kwargs)
-
-
-@pytest.fixture()
-def cftime_index(time_range_kwargs):
-    return xr.cftime_range(**time_range_kwargs)
+# Create a list of pairs of similar-length initial and resample frequencies
+# that cover:
+# - Resampling from shorter to longer frequencies
+# - Resampling from longer to shorter frequencies
+# - Resampling from one initial frequency to another.
+# These are used to test the cftime version of resample against pandas
+# with a standard calendar.
+FREQS = [
+    ('8003D', '4001D'),
+    ('8003D', '16006D'),
+    ('8003D', '21AS'),
+    ('6H', '3H'),
+    ('6H', '12H'),
+    ('6H', '400T'),
+    ('3D', 'D'),
+    ('3D', '6D'),
+    ('11D', 'MS'),
+    ('3MS', 'MS'),
+    ('3MS', '6MS'),
+    ('3MS', '85D'),
+    ('7M', '3M'),
+    ('7M', '14M'),
+    ('7M', '2QS-APR'),
+    ('43QS-AUG', '21QS-AUG'),
+    ('43QS-AUG', '86QS-AUG'),
+    ('43QS-AUG', '11A-JUN'),
+    ('11Q-JUN', '5Q-JUN'),
+    ('11Q-JUN', '22Q-JUN'),
+    ('11Q-JUN', '51MS'),
+    ('3AS-MAR', 'AS-MAR'),
+    ('3AS-MAR', '6AS-MAR'),
+    ('3AS-MAR', '14Q-FEB'),
+    ('7A-MAY', '3A-MAY'),
+    ('7A-MAY', '14A-MAY'),
+    ('7A-MAY', '85M')
+]
 
 
 def da(index):
@@ -38,65 +54,31 @@ def da(index):
                         coords=[index], dims=['time'])
 
 
-@pytest.mark.parametrize('freq', [
-    '700T', '8001T',
-    '12H', '8001H',
-    '8D', '8001D',
-    '2MS', '3MS',
-    '2QS-AUG', '3QS-SEP',
-    '3AS-MAR', '4A-MAY'])
-@pytest.mark.parametrize('closed', [None, 'right'])
-@pytest.mark.parametrize('label', [None, 'right'])
+@pytest.mark.parametrize('freqs', FREQS, ids=lambda x: '{}->{}'.format(*x))
+@pytest.mark.parametrize('closed', [None, 'left', 'right'])
+@pytest.mark.parametrize('label', [None, 'left', 'right'])
 @pytest.mark.parametrize('base', [24, 31])
-def test_resampler(freq, closed, label, base,
-                   datetime_index, cftime_index):
-    # Fairly extensive testing for standard/proleptic Gregorian calendar
-    # For any frequencies which are not greater-than-day and anchored
-    # at the end, the default values for closed and label are 'left'.
-    loffset = '12H'
-    try:
-        da_datetime = da(datetime_index).resample(
-            time=freq, closed=closed, label=label, base=base,
-            loffset=loffset).mean()
-    except ValueError:
-        with pytest.raises(ValueError):
-            da(cftime_index).resample(
-                time=freq, closed=closed, label=label, base=base,
-                loffset=loffset).mean()
-    else:
-        da_cftime = da(cftime_index).resample(time=freq, closed=closed,
-                                              label=label, base=base,
-                                              loffset=loffset).mean()
-        da_cftime['time'] = da_cftime.indexes['time'].to_datetimeindex()
-        xr.testing.assert_identical(da_cftime, da_datetime)
-
+def test_resample(freqs, closed, label, base):
+    initial_freq, resample_freq = freqs
+    start = '2000-01-01T12:07:01'
+    index_kwargs = dict(start=start, periods=5, freq=initial_freq)
+    datetime_index = pd.date_range(**index_kwargs)
+    cftime_index = xr.cftime_range(**index_kwargs)
 
-@pytest.mark.parametrize('freq', [
-    '2M', '3M',
-    '2Q-JUN', '3Q-JUL',
-    '3A-FEB', '4A-APR'])
-@pytest.mark.parametrize('closed', ['left', None])
-@pytest.mark.parametrize('label', ['left', None])
-@pytest.mark.parametrize('base', [17, 24])
-def test_resampler_end_super_day(freq, closed, label, base,
-                                 datetime_index, cftime_index):
-    # Fairly extensive testing for standard/proleptic Gregorian calendar.
-    # For greater-than-day frequencies anchored at the end, the default values
-    # for closed and label are 'right'.
     loffset = '12H'
     try:
         da_datetime = da(datetime_index).resample(
-            time=freq, closed=closed, label=label, base=base,
+            time=resample_freq, closed=closed, label=label, base=base,
             loffset=loffset).mean()
     except ValueError:
         with pytest.raises(ValueError):
             da(cftime_index).resample(
-                time=freq, closed=closed, label=label, base=base,
+                time=resample_freq, closed=closed, label=label, base=base,
                 loffset=loffset).mean()
     else:
-        da_cftime = da(cftime_index).resample(time=freq, closed=closed,
-                                              label=label, base=base,
-                                              loffset=loffset).mean()
+        da_cftime = da(cftime_index).resample(
+            time=resample_freq, closed=closed,
+            label=label, base=base, loffset=loffset).mean()
         da_cftime['time'] = da_cftime.indexes['time'].to_datetimeindex()
         xr.testing.assert_identical(da_cftime, da_datetime)
 
@@ -111,6 +93,7 @@ def test_closed_label_defaults(freq, expected):
     assert CFTimeGrouper(freq=freq).label == expected
 
 
+@pytest.mark.filterwarnings('ignore:Converting a CFTimeIndex')
 @pytest.mark.parametrize('calendar', ['gregorian', 'noleap', 'all_leap',
                                       '360_day', 'julian'])
 def test_calendars(calendar):
diff --git a/xarray/tests/test_combine.py b/xarray/tests/test_combine.py
index 0d03b6e0..6d0f4626 100644
--- a/xarray/tests/test_combine.py
+++ b/xarray/tests/test_combine.py
@@ -13,7 +13,7 @@ from xarray.core.combine import (
     _infer_tile_ids_from_nested_list, _new_tile_id)
 
 from . import (
-    InaccessibleArray, assert_array_equal, assert_combined_tile_ids_equal,
+    InaccessibleArray, assert_array_equal,
     assert_equal, assert_identical, raises_regex, requires_dask)
 from .test_dataset import create_test_data
 
@@ -418,6 +418,13 @@ class TestAutoCombine(object):
         assert_identical(expected, actual)
 
 
+def assert_combined_tile_ids_equal(dict1, dict2):
+    assert len(dict1) == len(dict2)
+    for k, v in dict1.items():
+        assert k in dict2.keys()
+        assert_equal(dict1[k], dict2[k])
+
+
 class TestTileIDsFromNestedList(object):
     def test_1d(self):
         ds = create_test_data
diff --git a/xarray/tests/test_dataset.py b/xarray/tests/test_dataset.py
index ab83d385..3ace80f5 100644
--- a/xarray/tests/test_dataset.py
+++ b/xarray/tests/test_dataset.py
@@ -2053,12 +2053,8 @@ class TestDataset(object):
         if python36_plus:
             with raises_regex(ValueError, 'both keyword and positional'):
                 original.expand_dims(OrderedDict((("d", 4),)), e=4)
-        else:
-            # In python 3.5, using dim_kwargs should raise a ValueError.
-            with raises_regex(ValueError, "dim_kwargs isn't"):
-                original.expand_dims(OrderedDict((("d", 4),)), e=4)
 
-    def test_expand_dims(self):
+    def test_expand_dims_int(self):
         original = Dataset({'x': ('a', np.random.randn(3)),
                             'y': (['b', 'a'], np.random.randn(4, 3))},
                            coords={'a': np.linspace(0, 1, 3),
@@ -2091,9 +2087,37 @@ class TestDataset(object):
         roundtripped = actual.squeeze('z')
         assert_identical(original, roundtripped)
 
+    def test_expand_dims_coords(self):
+        original = Dataset({'x': ('a', np.array([1, 2, 3]))})
+        expected = Dataset(
+            {'x': (('b', 'a'), np.array([[1, 2, 3], [1, 2, 3]]))},
+            coords={'b': [1, 2]},
+        )
+        actual = original.expand_dims(OrderedDict(b=[1, 2]))
+        assert_identical(expected, actual)
+        assert 'b' not in original._coord_names
+
+    def test_expand_dims_existing_scalar_coord(self):
+        original = Dataset({'x': 1}, {'a': 2})
+        expected = Dataset({'x': (('a',), [1])}, {'a': [2]})
+        actual = original.expand_dims('a')
+        assert_identical(expected, actual)
+
+    def test_isel_expand_dims_roundtrip(self):
+        original = Dataset({'x': (('a',), [1])}, {'a': [2]})
+        actual = original.isel(a=0).expand_dims('a')
+        assert_identical(actual, original)
+
+    def test_expand_dims_mixed_int_and_coords(self):
         # Test expanding one dimension to have size > 1 that doesn't have
         # coordinates, and also expanding another dimension to have size > 1
         # that DOES have coordinates.
+        original = Dataset({'x': ('a', np.random.randn(3)),
+                            'y': (['b', 'a'], np.random.randn(4, 3))},
+                           coords={'a': np.linspace(0, 1, 3),
+                                   'b': np.linspace(0, 1, 4),
+                                   'c': np.linspace(0, 1, 5)})
+
         actual = original.expand_dims(
             OrderedDict((("d", 4), ("e", ["l", "m", "n"]))))
 
@@ -2109,34 +2133,45 @@ class TestDataset(object):
                                            b=np.linspace(0, 1, 4),
                                            a=np.linspace(0, 1, 3)),
                                dims=['d', 'e', 'b', 'a']).drop('d')},
-            coords={'c': np.linspace(0, 1, 5)},
-            attrs={'key': 'entry'})
+            coords={'c': np.linspace(0, 1, 5)})
         assert_identical(actual, expected)
 
-        # Test with kwargs instead of passing dict to dim arg.
-
-        # TODO: only the code under the if-statement is needed when python 3.5
-        #   is no longer supported.
-        python36_plus = sys.version_info[0] == 3 and sys.version_info[1] > 5
-        if python36_plus:
-            other_way = original.expand_dims(e=["l", "m", "n"])
-            other_way_expected = Dataset(
-                {'x': xr.DataArray(original['x'].values * np.ones([3, 3]),
-                                   coords=dict(e=['l', 'm', 'n'],
-                                               a=np.linspace(0, 1, 3)),
-                                   dims=['e', 'a']),
-                 'y': xr.DataArray(original['y'].values * np.ones([3, 4, 3]),
-                                   coords=dict(e=['l', 'm', 'n'],
-                                               b=np.linspace(0, 1, 4),
-                                               a=np.linspace(0, 1, 3)),
-                                   dims=['e', 'b', 'a'])},
-                coords={'c': np.linspace(0, 1, 5)},
-                attrs={'key': 'entry'})
-            assert_identical(other_way_expected, other_way)
-        else:
-            # In python 3.5, using dim_kwargs should raise a ValueError.
-            with raises_regex(ValueError, "dim_kwargs isn't"):
-                original.expand_dims(e=["l", "m", "n"])
+    @pytest.mark.skipif(
+        sys.version_info[:2] > (3, 5),
+        reason="we only raise these errors for Python 3.5",
+    )
+    def test_expand_dims_kwargs_python35(self):
+        original = Dataset({'x': ('a', np.random.randn(3))})
+        with raises_regex(ValueError, "dim_kwargs isn't"):
+            original.expand_dims(e=["l", "m", "n"])
+        with raises_regex(TypeError, "must be an OrderedDict"):
+            original.expand_dims({'e': ["l", "m", "n"]})
+
+    @pytest.mark.skipif(
+        sys.version_info[:2] < (3, 6),
+        reason='keyword arguments are only ordered on Python 3.6+',
+    )
+    def test_expand_dims_kwargs_python36plus(self):
+        original = Dataset({'x': ('a', np.random.randn(3)),
+                            'y': (['b', 'a'], np.random.randn(4, 3))},
+                           coords={'a': np.linspace(0, 1, 3),
+                                   'b': np.linspace(0, 1, 4),
+                                   'c': np.linspace(0, 1, 5)},
+                           attrs={'key': 'entry'})
+        other_way = original.expand_dims(e=["l", "m", "n"])
+        other_way_expected = Dataset(
+            {'x': xr.DataArray(original['x'].values * np.ones([3, 3]),
+                               coords=dict(e=['l', 'm', 'n'],
+                                           a=np.linspace(0, 1, 3)),
+                               dims=['e', 'a']),
+             'y': xr.DataArray(original['y'].values * np.ones([3, 4, 3]),
+                               coords=dict(e=['l', 'm', 'n'],
+                                           b=np.linspace(0, 1, 4),
+                                           a=np.linspace(0, 1, 3)),
+                               dims=['e', 'b', 'a'])},
+            coords={'c': np.linspace(0, 1, 5)},
+            attrs={'key': 'entry'})
+        assert_identical(other_way_expected, other_way)
 
     def test_set_index(self):
         expected = create_test_multiindex()
diff --git a/xarray/tests/test_interp.py b/xarray/tests/test_interp.py
index 5596bfb3..8347d54b 100644
--- a/xarray/tests/test_interp.py
+++ b/xarray/tests/test_interp.py
@@ -291,7 +291,7 @@ def test_errors(use_dask):
     if use_dask:
         da = get_example_data(3)
     else:
-        da = get_example_data(1)
+        da = get_example_data(0)
 
     result = da.interp(x=[-1, 1, 3], kwargs={'fill_value': 0.0})
     assert not np.isnan(result.values).any()
diff --git a/xarray/tests/test_utils.py b/xarray/tests/test_utils.py
index e98ab5cd..caff78b4 100644
--- a/xarray/tests/test_utils.py
+++ b/xarray/tests/test_utils.py
@@ -5,15 +5,12 @@ import numpy as np
 import pandas as pd
 import pytest
 
-import xarray as xr
 from xarray.coding.cftimeindex import CFTimeIndex
 from xarray.core import duck_array_ops, utils
 from xarray.core.utils import either_dict_or_kwargs
-from xarray.testing import assert_identical
 
 from . import (
-    assert_array_equal, has_cftime, has_cftime_or_netCDF4, requires_cftime,
-    requires_dask)
+    assert_array_equal, has_cftime, has_cftime_or_netCDF4, requires_dask)
 from .test_coding_times import _all_cftime_date_types
 
 
@@ -178,19 +175,6 @@ class TestDictionaries(object):
         assert repr(utils.SortedKeysDict()) == \
             "SortedKeysDict({})"
 
-    def test_chain_map(self):
-        m = utils.ChainMap({'x': 0, 'y': 1}, {'x': -100, 'z': 2})
-        assert 'x' in m
-        assert 'y' in m
-        assert 'z' in m
-        assert m['x'] == 0
-        assert m['y'] == 1
-        assert m['z'] == 2
-        m['x'] = 100
-        assert m['x'] == 100
-        assert m.maps[0]['x'] == 100
-        assert set(m) == {'x', 'y', 'z'}
-
 
 def test_repr_object():
     obj = utils.ReprObject('foo')