Repository: neurostuff/PyMARE

Issue:
Datasets do not check input array shapes/sizes
To reproduce the problem, use the following code

```python
from pymare import core, estimators

y = [
    [2, 4, 6],  # estimates for first study's three datasets
    [3, 2, 1],  # estimates for second study's three datasets
]
v = [
    [100, 100, 100],  # estimate variance for first study's three datasets
    [8, 4, 2],  # estimate variance for second study's three datasets
]
X = [  # all "dataset"s must have the same regressors
    [5, 9],  # regressors for first study
    [2, 8],  # regressors for second study
    [5, 5],  # regressors for imaginary third study
]

dataset = core.Dataset(y=y, v=v, X=X, X_names=["X1", "X7"], add_intercept=False)
est = estimators.WeightedLeastSquares().fit_dataset(dataset)
results = est.summary()
```

Here's the traceback:

```Text
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-1-b48664cc0a71> in <module>
     16 
     17 dataset = core.Dataset(y=y, v=v, X=X, X_names=["X1", "X7"], add_intercept=False)
---> 18 est = estimators.WeightedLeastSquares().fit_dataset(dataset)
     19 results = est.summary()

~/Documents/tsalo/PyMARE/pymare/estimators/estimators.py in fit_dataset(self, dataset, *args, **kwargs)
     97 
     98         all_kwargs.update(kwargs)
---> 99         self.fit(*args, **all_kwargs)
    100         self.dataset_ = dataset
    101 

~/Documents/tsalo/PyMARE/pymare/estimators/estimators.py in fit(self, y, X, v)
    196             v = np.ones_like(y)
    197 
--> 198         beta, inv_cov = weighted_least_squares(y, v, X, self.tau2, return_cov=True)
    199         self.params_ = {"fe_params": beta, "tau2": self.tau2, "inv_cov": inv_cov}
    200         return self

~/Documents/tsalo/PyMARE/pymare/stats.py in weighted_least_squares(y, v, X, tau2, return_cov)
     33 
     34     # Einsum indices: k = studies, p = predictors, i = parallel iterates
---> 35     wX = np.einsum("kp,ki->ipk", X, w)
     36     cov = wX.dot(X)
     37 

<__array_function__ internals> in einsum(*args, **kwargs)

/opt/miniconda3/lib/python3.8/site-packages/numpy/core/einsumfunc.py in einsum(out, optimize, *operands, **kwargs)
   1357         if specified_out:
   1358             kwargs['out'] = out
-> 1359         return c_einsum(*operands, **kwargs)
   1360 
   1361     # Check the kwargs to avoid a more cryptic error later, without having to

ValueError: operands could not be broadcast together with remapped shapes [original->remapped]: (3,2)->(2,3) (2,3)->(3,newaxis,2) 
```

Note that the exception is raised when the Estimator is fitted, rather than when the Dataset is initialized, and it's the product of a mathematical operation, rather than an explicit check.

NOTE: Related to #37.

Discussion:
codecov-commenter: # [Codecov](https://codecov.io/gh/neurostuff/PyMARE/pull/106?src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=neurostuff) Report
> Merging [#106](https://codecov.io/gh/neurostuff/PyMARE/pull/106?src=pr&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=neurostuff) (04b7dbb) into [master](https://codecov.io/gh/neurostuff/PyMARE/commit/0a9984060b912c241633d74c9d5a91779a0a71f5?el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=neurostuff) (0a99840) will **decrease** coverage by `0.28%`.
> The diff coverage is `76.19%`.

```diff
@@            Coverage Diff             @@
##           master     #106      +/-   ##
==========================================
- Coverage   87.60%   87.31%   -0.29%     
==========================================
  Files          13       13              
  Lines         863      883      +20     
==========================================
+ Hits          756      771      +15     
- Misses        107      112       +5     
```


| [Impacted Files](https://codecov.io/gh/neurostuff/PyMARE/pull/106?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=neurostuff) | Coverage Δ | |
|---|---|---|
| [pymare/utils.py](https://codecov.io/gh/neurostuff/PyMARE/pull/106/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=neurostuff#diff-cHltYXJlL3V0aWxzLnB5) | `78.26% <70.58%> (-21.74%)` | :arrow_down: |
| [pymare/core.py](https://codecov.io/gh/neurostuff/PyMARE/pull/106/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=neurostuff#diff-cHltYXJlL2NvcmUucHk=) | `94.52% <100.00%> (+0.23%)` | :arrow_up: |

------

[Continue to review full report at Codecov](https://codecov.io/gh/neurostuff/PyMARE/pull/106?src=pr&el=continue&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=neurostuff).
> **Legend** - [Click here to learn more](https://docs.codecov.io/docs/codecov-delta?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=neurostuff)
> `Δ = absolute <relative> (impact)`, `ø = not affected`, `? = missing data`
> Powered by [Codecov](https://codecov.io/gh/neurostuff/PyMARE/pull/106?src=pr&el=footer&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=neurostuff). Last update [0a99840...04b7dbb](https://codecov.io/gh/neurostuff/PyMARE/pull/106?src=pr&el=lastupdated&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=neurostuff). Read the [comment docs](https://docs.codecov.io/docs/pull-request-comments?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=neurostuff).

tsalo: Oh BTW you should add yourself to the Zenodo file. I totally forgot about that.
JulioAPeraza: > I think we'll need to refactor _check_inputs_shape if we ever need to extend it to work with N-dimensional arrays

That's a great idea. I think that can be implemented by checking the shape looping through a list of axis given by the user:

```python
utils._check_inputs_shape(y, X, "y", "X", axis=[0])
utils._check_inputs_shape(y, n, "y", "n", axis=[0, 1])
utils._check_inputs_shape(X, np.array(X_names)[None, :], "X", "X_names", axis=[1])
```


Patch:
diff --git a/.zenodo.json b/.zenodo.json
index f524e90..00d9ae2 100644
--- a/.zenodo.json
+++ b/.zenodo.json
@@ -10,6 +10,11 @@
       "affiliation": "Florida International University",
       "orcid": "0000-0001-9813-3167"
     },
+    {
+      "name": "Peraza, Julio A.",
+      "affiliation": "Florida International University",
+      "orcid": "0000-0003-3816-5903"
+    },
     {
       "name": "Nichols, Thomas E.",
       "affiliation": "Big Data Institute, University of Oxford",
diff --git a/pymare/core.py b/pymare/core.py
index 74c5b48..7fc7b60 100644
--- a/pymare/core.py
+++ b/pymare/core.py
@@ -5,7 +5,7 @@ from functools import partial
 import numpy as np
 import pandas as pd
 
-from pymare.utils import _listify
+from pymare.utils import _check_inputs_shape, _listify
 
 from .estimators import (
     DerSimonianLaird,
@@ -94,6 +94,10 @@ class Dataset:
         self.X = X
         self.X_names = names
 
+        _check_inputs_shape(self.y, self.X, "y", "X", row=True)
+        _check_inputs_shape(self.y, self.v, "y", "v", row=True, column=True)
+        _check_inputs_shape(self.y, self.n, "y", "n", row=True, column=True)
+
     def _get_predictors(self, X, names, add_intercept):
         if X is None and not add_intercept:
             raise ValueError(
diff --git a/pymare/utils.py b/pymare/utils.py
index 1e11cb7..b159a48 100644
--- a/pymare/utils.py
+++ b/pymare/utils.py
@@ -19,3 +19,39 @@ def _listify(obj):
     This provides a simple way to accept flexible arguments.
     """
     return obj if isinstance(obj, (list, tuple, type(None), np.ndarray)) else [obj]
+
+
+def _check_inputs_shape(param1, param2, param1_name, param2_name, row=False, column=False):
+    """Check whether 'param1' and 'param2' have the same shape.
+
+    Parameters
+    ----------
+    param1 : array
+    param2 : array
+    param1_name : str
+    param2_name : str
+    row : bool, default to False.
+    column : bool, default to False.
+    """
+    if (param1 is not None) and (param2 is not None):
+        if row and not column:
+            shape1 = param1.shape[0]
+            shape2 = param2.shape[0]
+            message = "rows"
+        elif column and not row:
+            shape1 = param1.shape[1]
+            shape2 = param2.shape[1]
+            message = "columns"
+        elif row and column:
+            shape1 = param1.shape
+            shape2 = param2.shape
+            message = "rows and columns"
+        else:
+            raise ValueError("At least one of the two parameters (row or column) should be True.")
+
+        if shape1 != shape2:
+            raise ValueError(
+                f"{param1_name} and {param2_name} should have the same number of {message}. "
+                f"You provided {param1_name} with shape {param1.shape} and {param2_name} "
+                f"with shape {param2.shape}."
+            )


Test patch:
diff --git a/pymare/tests/test_utils.py b/pymare/tests/test_utils.py
index 43a81bf..474b8f8 100644
--- a/pymare/tests/test_utils.py
+++ b/pymare/tests/test_utils.py
@@ -1,6 +1,9 @@
 """Tests for pymare.utils."""
 import os.path as op
 
+import numpy as np
+import pytest
+
 from pymare import utils
 
 
@@ -8,3 +11,30 @@ def test_get_resource_path():
     """Test nimare.utils.get_resource_path."""
     print(utils.get_resource_path())
     assert op.isdir(utils.get_resource_path())
+
+
+def test_check_inputs_shape():
+    """Test nimare.utils._check_inputs_shape."""
+    n_rows = 5
+    n_columns = 4
+    n_pred = 3
+    y = np.random.randint(1, 100, size=(n_rows, n_columns))
+    v = np.random.randint(1, 100, size=(n_rows + 1, n_columns))
+    n = np.random.randint(1, 100, size=(n_rows, n_columns))
+    X = np.random.randint(1, 100, size=(n_rows, n_pred))
+    X_names = [f"X{x}" for x in range(n_pred)]
+
+    utils._check_inputs_shape(y, X, "y", "X", row=True)
+    utils._check_inputs_shape(y, n, "y", "n", row=True, column=True)
+    utils._check_inputs_shape(X, np.array(X_names)[None, :], "X", "X_names", column=True)
+
+    # Raise error if the number of rows and columns of v don't match y
+    with pytest.raises(ValueError):
+        utils._check_inputs_shape(y, v, "y", "v", row=True, column=True)
+
+    # Raise error if neither row or column is True
+    with pytest.raises(ValueError):
+        utils._check_inputs_shape(y, n, "y", "n")
+
+    # Dataset may be initialized with n or v as None
+    utils._check_inputs_shape(y, None, "y", "n", row=True, column=True)