CCA issues
There are actually three problems going on
1. The environmental labels are being dropped.  [Here](https://github.com/biocore/scikit-bio/blob/master/skbio/stats/ordination/_canonical_correspondence_analysis.py#L224), you can see that no column or row labels are being passed into biplot_scores
2. The column rows don't actually match up with what is given in vegan.  In the example from [SO](http://stackoverflow.com/questions/36208741/how-can-i-make-a-cca-triplot-using-scikit-bio-python), there are only 3 columns (i.e. 3 dimension). This probably should be consistent with the PC ids.  
3. Its also becoming apparent that there is a sore need for biplots.  CA, RDA and CCA are essentially unusable without this functionality present in scikit-bio.  Its not trivial to tune these plots by hand. 

diff --git a/CHANGELOG.md b/CHANGELOG.md
index a86e8d592..37f1879e8 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -18,6 +18,7 @@
 * `TabularMSA.append` and `TabularMSA.extend` now require one of `minter`, `index`, or `reset_index` to be provided when incorporating new sequences into an MSA. Previous behavior was to auto-increment the index labels if `minter` and `index` weren't provided and the MSA had a default integer index, otherwise error. Use `reset_index=True` to obtain the previous behavior in a more explicit way.
 
 ### Bug fixes
+* Fixed row and column names to `biplot_scores` in the `OrdinationResults` object from `skbio.stats.ordination`. This fix affect the `cca` and `rda` methods. ([#1322](https://github.com/biocore/scikit-bio/issues/1322))
 * Fixed bug when using `skbio.io.format.stockholm` reader on file with multi-line tree with no id. Previously this raised an `AttributeError`, now it correctly handles this type of tree. ([#1334](https://github.com/biocore/scikit-bio/issues/1334))
 * Fixed bug when reading Stockholm files with GF or GS features split over multiple lines. Previously, the feature text was simply concatenated because it was assumed to have trailing whitespace. There are examples of Stockholm files with and without trailing whitespace for multi-line features, so the `skbio.io.format.stockholm` reader now adds a single space when concatenating feature text without trailing whitespace to avoid joining words together. Multi-line trees stored as GF metadata are concatenated as they appear in the file; a space is not added when concatenating. ([#1328](https://github.com/biocore/scikit-bio/issues/1328))
 * Fixed bug when using `Sequence.iter_kmers` on empty `Sequence` object. Previously this raised a `ValueError`, now it returns
diff --git a/skbio/io/format/tests/test_ordination.py b/skbio/io/format/tests/test_ordination.py
index ad8961676..2c37a7d32 100644
--- a/skbio/io/format/tests/test_ordination.py
+++ b/skbio/io/format/tests/test_ordination.py
@@ -194,7 +194,7 @@ class OrdinationResultsReaderWriterTests(OrdinationTestData):
                 obs = _ordination_to_ordination_results(fp)
                 assert_ordination_results_equal(
                     obs, obj, ignore_method_names=True,
-                    ignore_axis_labels=True, ignore_biplot_scores_labels=True)
+                    ignore_axis_labels=True)
 
     def test_read_invalid_files(self):
         for invalid_fp, error_msg_regexp, _ in self.invalid_fps:
diff --git a/skbio/stats/ordination/_canonical_correspondence_analysis.py b/skbio/stats/ordination/_canonical_correspondence_analysis.py
index 5ba436dde..18b4265f1 100644
--- a/skbio/stats/ordination/_canonical_correspondence_analysis.py
+++ b/skbio/stats/ordination/_canonical_correspondence_analysis.py
@@ -221,7 +221,10 @@ def cca(y, x, scaling=1):
                            columns=pc_ids, index=sample_ids)
     features = pd.DataFrame(features_scores,
                             columns=pc_ids, index=feature_ids)
-    biplot_scores = pd.DataFrame(biplot_scores)
+
+    biplot_scores = pd.DataFrame(biplot_scores,
+                                 index=x.columns,
+                                 columns=pc_ids[:biplot_scores.shape[1]])
     sample_constraints = pd.DataFrame(sample_constraints,
                                       index=sample_ids, columns=pc_ids)
 
diff --git a/skbio/stats/ordination/_redundancy_analysis.py b/skbio/stats/ordination/_redundancy_analysis.py
index e5ef21339..252964eb4 100644
--- a/skbio/stats/ordination/_redundancy_analysis.py
+++ b/skbio/stats/ordination/_redundancy_analysis.py
@@ -201,7 +201,10 @@ def rda(y, x, scale_Y=False, scaling=1):
     # can see that there's an arrow for each of the 4
     # environmental variables (depth, coral, sand, other) even if
     # other = not(coral or sand)
-    biplot_scores = pd.DataFrame(corr(X, u))
+    biplot_scores = corr(X, u)
+    biplot_scores = pd.DataFrame(biplot_scores,
+                                 index=x.columns,
+                                 columns=pc_ids[:biplot_scores.shape[1]])
     # The "Correlations of environmental variables with sample
     # scores" from table 11.4 are quite similar to vegan's biplot
     # scores, but they're computed like this:
diff --git a/skbio/stats/ordination/tests/test_canonical_correspondence_analysis.py b/skbio/stats/ordination/tests/test_canonical_correspondence_analysis.py
index 78778f2bb..b21944671 100644
--- a/skbio/stats/ordination/tests/test_canonical_correspondence_analysis.py
+++ b/skbio/stats/ordination/tests/test_canonical_correspondence_analysis.py
@@ -59,7 +59,7 @@ class TestCCAResults1(TestCase):
                            'Sample4', 'Sample5', 'Sample6', 'Sample7',
                            'Sample8', 'Sample9']
         self.env_ids = ['Constraint0', 'Constraint1',
-                        'Constraint2', 'Constraint3']
+                        'Constraint2']
         self.pc_ids = ['CCA1', 'CCA2', 'CCA3', 'CCA4', 'CCA5', 'CCA6', 'CCA7',
                        'CCA8', 'CCA9']
         self.Y = pd.DataFrame(
@@ -67,10 +67,10 @@ class TestCCAResults1(TestCase):
             columns=self.feature_ids,
             index=self.sample_ids)
         self.X = pd.DataFrame(
-            np.loadtxt(get_data_path('example3_X')),
+            np.loadtxt(get_data_path('example3_X'))[:, :-1],
             columns=self.env_ids,
             index=self.sample_ids
-            ).ix[:, :-1]
+            )
 
     def test_scaling1(self):
         scores = cca(self.Y, self.X, scaling=1)
@@ -93,10 +93,12 @@ class TestCCAResults1(TestCase):
                 'example3_sample_constraints_scaling1')),
             index=self.sample_ids,
             columns=self.pc_ids)
-
-        biplot_scores = pd.DataFrame(
-            np.loadtxt(get_data_path(
-                'example3_biplot_scaling1')))
+        mat = np.loadtxt(get_data_path(
+            'example3_biplot_scaling1'))
+        cropped_pcs = self.pc_ids[:mat.shape[1]]
+        biplot_scores = pd.DataFrame(mat,
+                                     index=self.env_ids,
+                                     columns=cropped_pcs)
 
         proportion_explained = pd.Series([0.466911, 0.238327, 0.100548,
                                           0.104937, 0.044805, 0.029747,
@@ -116,7 +118,6 @@ class TestCCAResults1(TestCase):
             eigvals=eigvals)
 
         assert_ordination_results_equal(scores, exp,
-                                        ignore_biplot_scores_labels=True,
                                         decimal=6)
 
     def test_scaling2(self):
@@ -141,9 +142,13 @@ class TestCCAResults1(TestCase):
             index=self.sample_ids,
             columns=self.pc_ids)
 
-        biplot_scores = pd.DataFrame(
-            np.loadtxt(get_data_path(
-                'example3_biplot_scaling2')))
+        mat = np.loadtxt(get_data_path(
+            'example3_biplot_scaling2'))
+
+        cropped_pc_ids = self.pc_ids[:mat.shape[1]]
+        biplot_scores = pd.DataFrame(mat,
+                                     index=self.env_ids,
+                                     columns=cropped_pc_ids)
 
         proportion_explained = pd.Series([0.466911, 0.238327, 0.100548,
                                           0.104937, 0.044805, 0.029747,
@@ -163,7 +168,6 @@ class TestCCAResults1(TestCase):
             eigvals=eigvals)
 
         assert_ordination_results_equal(scores, exp,
-                                        ignore_biplot_scores_labels=True,
                                         decimal=6)
 
 
diff --git a/skbio/stats/ordination/tests/test_redundancy_analysis.py b/skbio/stats/ordination/tests/test_redundancy_analysis.py
index 13c8a0e7a..c997a27f3 100644
--- a/skbio/stats/ordination/tests/test_redundancy_analysis.py
+++ b/skbio/stats/ordination/tests/test_redundancy_analysis.py
@@ -38,7 +38,7 @@ class TestRDAResults(TestCase):
                            'Site5', 'Site6', 'Site7', 'Site8', 'Site9']
         self.feature_ids = ['Species0', 'Species1', 'Species2', 'Species3',
                             'Species4', 'Species5']
-        self.env_ids = map(str, range(4))
+        self.env_ids = list(map(str, range(4)))
         self.pc_ids = ['RDA1', 'RDA2', 'RDA3', 'RDA4', 'RDA5', 'RDA6', 'RDA7']
 
         self.Y = pd.DataFrame(
@@ -53,9 +53,6 @@ class TestRDAResults(TestCase):
 
         scores = rda(self.Y, self.X, scaling=1)
 
-        biplot_scores = pd.DataFrame(np.loadtxt(
-            get_data_path('example2_biplot_scaling1')))
-
         sample_constraints = pd.DataFrame(np.loadtxt(
             get_data_path('example2_sample_constraints_scaling1')))
 
@@ -77,10 +74,12 @@ class TestRDAResults(TestCase):
                 'example2_sample_constraints_scaling1')),
             index=self.sample_ids,
             columns=self.pc_ids)
-
-        biplot_scores = pd.DataFrame(
-            np.loadtxt(get_data_path(
-                'example2_biplot_scaling1')))
+        mat = np.loadtxt(get_data_path(
+            'example2_biplot_scaling1'))
+        cropped_pc_ids = self.pc_ids[:mat.shape[1]]
+        biplot_scores = pd.DataFrame(mat,
+                                     index=self.env_ids,
+                                     columns=cropped_pc_ids)
 
         # These are wrong. See issue #1002
         proportion_explained = pd.Series([0.44275783, 0.25614586,
@@ -104,15 +103,16 @@ class TestRDAResults(TestCase):
 
         assert_ordination_results_equal(scores, exp,
                                         ignore_directionality=True,
-                                        ignore_biplot_scores_labels=True,
                                         decimal=6)
 
     def test_scaling2(self):
 
         scores = rda(self.Y, self.X, scaling=2)
-
-        biplot_scores = pd.DataFrame(np.loadtxt(
-            get_data_path('example2_biplot_scaling2')))
+        mat = np.loadtxt(get_data_path('example2_biplot_scaling2'))
+        cropped_pc_ids = self.pc_ids[:mat.shape[1]]
+        biplot_scores = pd.DataFrame(mat,
+                                     index=self.env_ids,
+                                     columns=cropped_pc_ids)
 
         sample_constraints = pd.DataFrame(np.loadtxt(
             get_data_path('example2_sample_constraints_scaling2')))
@@ -136,9 +136,12 @@ class TestRDAResults(TestCase):
             index=self.sample_ids,
             columns=self.pc_ids)
 
-        biplot_scores = pd.DataFrame(
-            np.loadtxt(get_data_path(
-                'example2_biplot_scaling2')))
+        mat = np.loadtxt(get_data_path(
+            'example2_biplot_scaling2'))
+        cropped_pc_ids = self.pc_ids[:mat.shape[1]]
+        biplot_scores = pd.DataFrame(mat,
+                                     index=self.env_ids,
+                                     columns=cropped_pc_ids)
 
         # These are wrong. See issue #1002
         proportion_explained = pd.Series([0.44275783, 0.25614586,
@@ -162,7 +165,6 @@ class TestRDAResults(TestCase):
 
         assert_ordination_results_equal(scores, exp,
                                         ignore_directionality=True,
-                                        ignore_biplot_scores_labels=True,
                                         decimal=6)
 
 
diff --git a/skbio/util/_testing.py b/skbio/util/_testing.py
index 708946a76..04465ddcd 100644
--- a/skbio/util/_testing.py
+++ b/skbio/util/_testing.py
@@ -150,7 +150,6 @@ def get_data_path(fn, subfolder='data'):
 @experimental(as_of="0.4.0")
 def assert_ordination_results_equal(left, right, ignore_method_names=False,
                                     ignore_axis_labels=False,
-                                    ignore_biplot_scores_labels=False,
                                     ignore_directionality=False,
                                     decimal=7):
     """Assert that ordination results objects are equal.
@@ -166,8 +165,6 @@ def assert_ordination_results_equal(left, right, ignore_method_names=False,
         Ignore differences in `short_method_name` and `long_method_name`.
     ignore_axis_labels : bool, optional
         Ignore differences in axis labels (i.e., column labels).
-    ignore_biplot_scores_labels : bool, optional
-        Ignore differences in `biplot_scores` row and column labels.
     ignore_directionality : bool, optional
         Ignore differences in directionality (i.e., differences in signs) for
         attributes `samples`, `features` and `biplot_scores`.
@@ -195,8 +192,7 @@ def assert_ordination_results_equal(left, right, ignore_method_names=False,
                         decimal=decimal)
 
     _assert_frame_equal(left.biplot_scores, right.biplot_scores,
-                        ignore_biplot_scores_labels,
-                        ignore_biplot_scores_labels,
+                        ignore_columns=ignore_axis_labels,
                         ignore_directionality=ignore_directionality,
                         decimal=decimal)
 
['skbio/stats/ordination/tests/test_redundancy_analysis.py', 'skbio/util/_testing.py', 'skbio/stats/ordination/_redundancy_analysis.py', 'skbio/stats/ordination/_canonical_correspondence_analysis.py', 'CHANGELOG.md', 'skbio/io/format/tests/test_ordination.py', 'skbio/stats/ordination/tests/test_canonical_correspondence_analysis.py']