diff --git a/.github/workflows/examples.yml b/.github/workflows/examples.yml
new file mode 100644
index 00000000..072b4b5a
--- /dev/null
+++ b/.github/workflows/examples.yml
@@ -0,0 +1,33 @@
+name: Examples
+
+env:
+  PROJECT_NAME: pfhedge
+  PYTHON_VERSION: '3.9'
+
+on:
+  workflow_dispatch:
+
+jobs:
+
+  examples:
+
+    name: Examples
+
+    runs-on: ubuntu-latest
+
+    steps:
+
+    - uses: actions/checkout@v2
+
+    - uses: actions/setup-python@v2
+      with:
+        python-version: ${{ env.PYTHON_VERSION }}
+
+    - name: Install
+      run: |
+        pip install poetry
+        poetry install
+        poetry run pip install -r examples/requirements.txt
+
+    - name: Run examples
+      run: cd examples && find . -name '*.py' -exec poetry run python {} \;
diff --git a/.gitignore b/.gitignore
index c3fc748f..a05d6b61 100644
--- a/.gitignore
+++ b/.gitignore
@@ -1,7 +1,7 @@
 poetry.lock
-docs/source/generated/
-examples/output/
-examples/output/.gitkeep
+docs/source/generated
+examples/output/*
+!examples/output/.gitkeep
 
 # Byte-compiled / optimized / DLL files
 __pycache__/
diff --git a/tests/nn/bs/__init__.py b/examples/output/.gitkeep
similarity index 100%
rename from tests/nn/bs/__init__.py
rename to examples/output/.gitkeep
diff --git a/pfhedge/_utils/parse.py b/pfhedge/_utils/parse.py
index 7c5cd7a4..1ed5c3b7 100644
--- a/pfhedge/_utils/parse.py
+++ b/pfhedge/_utils/parse.py
@@ -42,3 +42,12 @@ def parse_volatility(
         return variance.clamp(min=0.0).sqrt()
     else:
         raise ValueError("Insufficient parameters to parse volatility")
+
+
+def parse_time_to_maturity(
+    *, time_to_maturity: Optional[Tensor] = None, **kwargs
+) -> Tensor:
+    if time_to_maturity is not None:
+        return time_to_maturity
+    else:
+        raise ValueError("Insufficient parameters to parse time_to_maturity")
diff --git a/pfhedge/autogreek.py b/pfhedge/autogreek.py
index a23f02a6..2e0d5d88 100644
--- a/pfhedge/autogreek.py
+++ b/pfhedge/autogreek.py
@@ -5,6 +5,7 @@
 from torch import Tensor
 
 from ._utils.parse import parse_spot
+from ._utils.parse import parse_time_to_maturity
 from ._utils.parse import parse_volatility
 
 
@@ -231,3 +232,62 @@ def vega(
         grad_outputs=torch.ones_like(price),
         create_graph=create_graph,
     )[0]
+
+
+def theta(
+    pricer: Callable[..., Tensor], *, create_graph: bool = False, **params
+) -> Tensor:
+    """Computes and returns theta of a derivative using automatic differentiation.
+
+    Theta is a differentiation of a derivative price with respect to
+    a variance of underlying instrument.
+
+    Note:
+        The keyword argument ``**params`` should contain at least one of the
+        following parameters:
+
+        - ``time_to_maturity``
+
+    Args:
+        pricer (callable): Pricing formula of a derivative.
+        create_graph (bool, default=False): If ``True``,
+            graph of the derivative will be constructed,
+            allowing to compute higher order derivative products.
+        **params: Parameters passed to ``pricer``.
+
+    Returns:
+        torch.Tensor
+
+    Examples:
+
+        Theta of a European option can be evaluated as follows.
+
+        >>> import pfhedge.autogreek as autogreek
+        >>> from pfhedge.nn import BSEuropeanOption
+        >>>
+        >>> pricer = BSEuropeanOption().price
+        >>> autogreek.theta(
+        ...     pricer,
+        ...     log_moneyness=torch.zeros(3),
+        ...     time_to_maturity=torch.tensor([0.1, 0.2, 0.3]),
+        ...     volatility=torch.tensor([0.20, 0.20, 0.20]),
+        ... )
+        tensor([-0.1261, -0.0891, -0.0727])
+    """
+    time_to_maturity = parse_time_to_maturity(**params).requires_grad_()
+    params["time_to_maturity"] = time_to_maturity
+
+    # Delete parameters that are not in the signature of pricer to avoid
+    # TypeError: <pricer> got an unexpected keyword argument '<parameter>'
+    for parameter in list(params.keys()):
+        if parameter not in signature(pricer).parameters.keys():
+            del params[parameter]
+
+    price = pricer(**params)
+    # Note: usually theta is calculated reversely (\partial{S}/\partial{T} = \partial{S}/\partial{-time_to_maturity})
+    return -torch.autograd.grad(
+        price,
+        inputs=time_to_maturity,
+        grad_outputs=torch.ones_like(price),
+        create_graph=create_graph,
+    )[0]
diff --git a/pfhedge/nn/modules/bs/_base.py b/pfhedge/nn/modules/bs/_base.py
index b8e8924f..c9e15c3e 100644
--- a/pfhedge/nn/modules/bs/_base.py
+++ b/pfhedge/nn/modules/bs/_base.py
@@ -69,6 +69,15 @@ def vega(self, **kwargs) -> Tensor:
         """
         return autogreek.vega(self.price, **kwargs)
 
+    @no_type_check
+    def theta(self, **kwargs) -> Tensor:
+        """Returns delta of the derivative.
+
+        Returns:
+            torch.Tensor
+        """
+        return autogreek.theta(self.price, **kwargs)
+
     def inputs(self) -> List[str]:
         """Returns the names of input features.
 
diff --git a/pfhedge/nn/modules/bs/american_binary.py b/pfhedge/nn/modules/bs/american_binary.py
index 0506d85b..3919a9a4 100644
--- a/pfhedge/nn/modules/bs/american_binary.py
+++ b/pfhedge/nn/modules/bs/american_binary.py
@@ -203,6 +203,79 @@ def gamma(
             volatility=volatility,
         )
 
+    @torch.enable_grad()
+    def vega(
+        self,
+        log_moneyness: Tensor,
+        max_log_moneyness: Tensor,
+        time_to_maturity: Tensor,
+        volatility: Tensor,
+    ) -> Tensor:
+        """Returns vega of the derivative.
+
+        Args:
+            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor): Time to expiry of the option.
+            volatility (torch.Tensor): Volatility of the underlying asset.
+
+        Shape:
+            - log_moneyness: :math:`(N, *)` where
+              :math:`*` means any number of additional dimensions.
+            - max_log_moneyness: :math:`(N, *)`
+            - time_to_maturity: :math:`(N, *)`
+            - volatility: :math:`(N, *)`
+            - output: :math:`(N, *)`
+
+        Returns:
+            torch.Tensor
+        """
+        return super().vega(
+            strike=self.strike,
+            log_moneyness=log_moneyness,
+            max_log_moneyness=max_log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+        )
+
+    @torch.enable_grad()
+    def theta(
+        self,
+        log_moneyness: Tensor,
+        max_log_moneyness: Tensor,
+        time_to_maturity: Tensor,
+        volatility: Tensor,
+    ) -> Tensor:
+        """Returns theta of the derivative.
+
+        Args:
+            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor): Time to expiry of the option.
+            volatility (torch.Tensor): Volatility of the underlying asset.
+
+        Shape:
+            - log_moneyness: :math:`(N, *)` where
+              :math:`*` means any number of additional dimensions.
+            - max_log_moneyness: :math:`(N, *)`
+            - time_to_maturity: :math:`(N, *)`
+            - volatility: :math:`(N, *)`
+            - output: :math:`(N, *)`
+
+        Note:
+            Risk-free rate is set to zero.
+
+        Returns:
+            torch.Tensor
+        """
+        return super().theta(
+            strike=self.strike,
+            log_moneyness=log_moneyness,
+            max_log_moneyness=max_log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+        )
+
     def implied_volatility(
         self,
         log_moneyness: Tensor,
diff --git a/pfhedge/nn/modules/bs/european.py b/pfhedge/nn/modules/bs/european.py
index e4ddd159..835887bf 100644
--- a/pfhedge/nn/modules/bs/european.py
+++ b/pfhedge/nn/modules/bs/european.py
@@ -137,10 +137,6 @@ def gamma(
         Returns:
             torch.Tensor
         """
-        if not self.call:
-            raise ValueError(
-                f"{self.__class__.__name__} for a put option is not yet supported."
-            )
 
         s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
         price = self.strike * s.exp()
@@ -148,6 +144,60 @@ def gamma(
 
         return gamma
 
+    def vega(
+        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+    ) -> Tensor:
+        """Returns vega of the derivative.
+
+        Args:
+            log_moneyness: (torch.Tensor): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor): Time to expiry of the option.
+            volatility (torch.Tensor): Volatility of the underlying asset.
+
+        Shape:
+            - log_moneyness: :math:`(N, *)` where
+              :math:`*` means any number of additional dimensions.
+            - time_to_maturity: :math:`(N, *)`
+            - volatility: :math:`(N, *)`
+            - output: :math:`(N, *)`
+
+        Returns:
+            torch.Tensor
+        """
+        s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
+        price = self.strike * s.exp()
+        vega = npdf(d1(s, t, v)) * price * t.sqrt()
+
+        return vega
+
+    def theta(
+        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+    ) -> Tensor:
+        """Returns theta of the derivative.
+
+        Args:
+            log_moneyness: (torch.Tensor): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor): Time to expiry of the option.
+            volatility (torch.Tensor): Volatility of the underlying asset.
+
+        Shape:
+            - log_moneyness: :math:`(N, *)` where
+              :math:`*` means any number of additional dimensions.
+            - time_to_maturity: :math:`(N, *)`
+            - volatility: :math:`(N, *)`
+            - output: :math:`(N, *)`
+
+        Note:
+            Risk-free rate is set to zero.
+
+        Returns:
+            torch.Tensor
+        """
+        s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
+        price = self.strike * s.exp()
+        theta = -npdf(d1(s, t, v)) * price * v / (2 * t.sqrt())
+        return theta
+
     def price(
         self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
     ) -> Tensor:
diff --git a/pfhedge/nn/modules/bs/european_binary.py b/pfhedge/nn/modules/bs/european_binary.py
index 92b7d9ec..ef0d1408 100644
--- a/pfhedge/nn/modules/bs/european_binary.py
+++ b/pfhedge/nn/modules/bs/european_binary.py
@@ -177,6 +177,63 @@ def gamma(
             volatility=volatility,
         )
 
+    def vega(
+        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+    ) -> Tensor:
+        """Returns vega of the derivative.
+
+        Args:
+            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor): Time to expiry of the option.
+            volatility (torch.Tensor): Volatility of the underlying asset.
+
+        Shape:
+            - log_moneyness: :math:`(N, *)`
+            - time_to_maturity: :math:`(N, *)`
+            - volatility: :math:`(N, *)`
+            - output: :math:`(N, *)`
+
+        Returns:
+            torch.Tensor
+        """
+        # TODO: Directly compute theta.
+        return super().vega(
+            strike=self.strike,
+            log_moneyness=log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+        )
+
+    def theta(
+        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+    ) -> Tensor:
+        """Returns theta of the derivative.
+
+        Args:
+            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor): Time to expiry of the option.
+            volatility (torch.Tensor): Volatility of the underlying asset.
+
+        Shape:
+            - log_moneyness: :math:`(N, *)`
+            - time_to_maturity: :math:`(N, *)`
+            - volatility: :math:`(N, *)`
+            - output: :math:`(N, *)`
+
+        Note:
+            Risk-free rate is set to zero.
+
+        Returns:
+            torch.Tensor
+        """
+        # TODO: Directly compute theta.
+        return super().theta(
+            strike=self.strike,
+            log_moneyness=log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+        )
+
     def implied_volatility(
         self,
         log_moneyness: Tensor,
diff --git a/pfhedge/nn/modules/bs/lookback.py b/pfhedge/nn/modules/bs/lookback.py
index aa3cc176..72d83551 100644
--- a/pfhedge/nn/modules/bs/lookback.py
+++ b/pfhedge/nn/modules/bs/lookback.py
@@ -241,6 +241,81 @@ def gamma(
             volatility=volatility,
         )
 
+    @torch.enable_grad()
+    def vega(
+        self,
+        log_moneyness: Tensor,
+        max_log_moneyness: Tensor,
+        time_to_maturity: Tensor,
+        volatility: Tensor,
+    ) -> Tensor:
+        """Returns vega of the derivative.
+
+        Args:
+            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor): Time to expiry of the option.
+            volatility (torch.Tensor):
+                Volatility of the underlying asset.
+
+        Shape:
+            - log_moneyness: :math:`(N, *)` where
+              :math:`*` means any number of additional dimensions.
+            - max_log_moneyness: :math:`(N, *)`
+            - time_to_maturity: :math:`(N, *)`
+            - volatility: :math:`(N, *)`
+            - output: :math:`(N, *)`
+
+        Returns:
+            torch.Tensor
+        """
+        return super().vega(
+            strike=self.strike,
+            log_moneyness=log_moneyness,
+            max_log_moneyness=max_log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+        )
+
+    @torch.enable_grad()
+    def theta(
+        self,
+        log_moneyness: Tensor,
+        max_log_moneyness: Tensor,
+        time_to_maturity: Tensor,
+        volatility: Tensor,
+    ) -> Tensor:
+        """Returns theta of the derivative.
+
+        Args:
+            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor): Time to expiry of the option.
+            volatility (torch.Tensor):
+                Volatility of the underlying asset.
+
+        Shape:
+            - log_moneyness: :math:`(N, *)` where
+              :math:`*` means any number of additional dimensions.
+            - max_log_moneyness: :math:`(N, *)`
+            - time_to_maturity: :math:`(N, *)`
+            - volatility: :math:`(N, *)`
+            - output: :math:`(N, *)`
+
+        Note:
+            Risk-free rate is set to zero.
+
+        Returns:
+            torch.Tensor
+        """
+        return super().theta(
+            strike=self.strike,
+            log_moneyness=log_moneyness,
+            max_log_moneyness=max_log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+        )
+
     def implied_volatility(
         self,
         log_moneyness: Tensor,
diff --git a/pyproject.toml b/pyproject.toml
index 745fa144..b4e657a3 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -1,6 +1,6 @@
 [tool.poetry]
 name = "pfhedge"
-version = "0.14.2"
+version = "0.15.0"
 description = "Deep Hedging in PyTorch"
 authors = ["Shota Imaki <shota.imaki.0801@gmail.com>"]
 license = "MIT"
@@ -8,7 +8,7 @@ repository = "https://github.com/pfnet-research/pfhedge"
 
 [tool.poetry.dependencies]
 python = "^3.7.12"
-torch = "^1.10.0"
+torch = "^1.9.0"
 tqdm = "^4.62.3"
 
 [tool.poetry.dev-dependencies]
diff --git a/tests/_utils/test_parse.py b/tests/_utils/test_parse.py
index fe08e3f5..5f136538 100644
--- a/tests/_utils/test_parse.py
+++ b/tests/_utils/test_parse.py
@@ -3,6 +3,7 @@
 from torch.testing import assert_close
 
 from pfhedge._utils.parse import parse_spot
+from pfhedge._utils.parse import parse_time_to_maturity
 from pfhedge._utils.parse import parse_volatility
 
 
@@ -46,3 +47,17 @@ def test_parse_volatility():
         _ = parse_volatility()
     with pytest.raises(ValueError):
         _ = parse_volatility(spot=volatility)
+
+
+def test_parse_time_to_maturity():
+    torch.manual_seed(42)
+
+    time_to_maturity = torch.randn(10).exp()
+
+    result = parse_time_to_maturity(time_to_maturity=time_to_maturity)
+    assert_close(result, time_to_maturity)
+
+    with pytest.raises(ValueError):
+        _ = parse_time_to_maturity()
+    with pytest.raises(ValueError):
+        _ = parse_time_to_maturity(spot=time_to_maturity)
diff --git a/tests/nn/modules/bs/__init__.py b/tests/nn/modules/bs/__init__.py
new file mode 100644
index 00000000..e69de29b
diff --git a/tests/nn/bs/_base.py b/tests/nn/modules/bs/_base.py
similarity index 89%
rename from tests/nn/bs/_base.py
rename to tests/nn/modules/bs/_base.py
index 60a95486..cfaefbfe 100644
--- a/tests/nn/bs/_base.py
+++ b/tests/nn/modules/bs/_base.py
@@ -28,6 +28,12 @@ def assert_shape_delta(self, module):
     def assert_shape_gamma(self, module):
         self._assert_shape(module, "gamma")
 
+    def assert_shape_vega(self, module):
+        self._assert_shape(module, "vega")
+
+    def assert_shape_theta(self, module):
+        self._assert_shape(module, "theta")
+
     def assert_shape_price(self, module):
         self._assert_shape(module, "price")
 
diff --git a/tests/nn/bs/_utils.py b/tests/nn/modules/bs/_utils.py
similarity index 61%
rename from tests/nn/bs/_utils.py
rename to tests/nn/modules/bs/_utils.py
index 8859059b..07b10f71 100644
--- a/tests/nn/bs/_utils.py
+++ b/tests/nn/modules/bs/_utils.py
@@ -9,5 +9,13 @@ def compute_gamma(module, input: Tensor) -> Tensor:
     return module.gamma(*(input[..., i] for i in range(input.size(-1))))
 
 
+def compute_vega(module, input: Tensor) -> Tensor:
+    return module.vega(*(input[..., i] for i in range(input.size(-1))))
+
+
+def compute_theta(module, input: Tensor) -> Tensor:
+    return module.theta(*(input[..., i] for i in range(input.size(-1))))
+
+
 def compute_price(module, input: Tensor) -> Tensor:
     return module.price(*(input[..., i] for i in range(input.size(-1))))
diff --git a/tests/nn/bs/test_american_binary.py b/tests/nn/modules/bs/test_american_binary.py
similarity index 78%
rename from tests/nn/bs/test_american_binary.py
rename to tests/nn/modules/bs/test_american_binary.py
index 5af16451..941fc307 100644
--- a/tests/nn/bs/test_american_binary.py
+++ b/tests/nn/modules/bs/test_american_binary.py
@@ -14,6 +14,8 @@
 from ._utils import compute_delta
 from ._utils import compute_gamma
 from ._utils import compute_price
+from ._utils import compute_theta
+from ._utils import compute_vega
 
 
 class TestBSAmericanBinaryOption(_TestBSModule):
@@ -83,6 +85,52 @@ def test_check_gamma(self):
         expect = torch.tensor([0.0])
         assert_close(result, expect)
 
+    def test_check_vega(self):
+        m = BSAmericanBinaryOption()
+
+        # vega = 0 for max --> +0
+        result = compute_vega(m, torch.tensor([[-10.0, -10.0, 0.1, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # vega = 0 for max > 0
+        result = compute_vega(m, torch.tensor([[0.0, 0.01, 0.1, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # vega = 0 for time --> +0
+        result = compute_vega(m, torch.tensor([[-0.01, -0.01, 1e-10, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # vega = 0 for volatility --> +0
+        result = compute_vega(m, torch.tensor([[-0.01, -0.01, 0.1, 1e-10]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+    def test_check_theta(self):
+        m = BSAmericanBinaryOption()
+
+        # vega = 0 for max --> +0
+        result = compute_theta(m, torch.tensor([[-10.0, -10.0, 0.1, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # vega = 0 for max > 0
+        result = compute_theta(m, torch.tensor([[0.0, 0.01, 0.1, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # vega = 0 for time --> +0
+        result = compute_theta(m, torch.tensor([[-0.01, -0.01, 1e-10, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # vega = 0 for volatility --> +0
+        result = compute_theta(m, torch.tensor([[-0.01, -0.01, 0.1, 1e-10]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
     def test_check_price(self):
         m = BSAmericanBinaryOption()
 
@@ -173,5 +221,7 @@ def test_shape(self):
         m = BSAmericanBinaryOption()
         self.assert_shape_delta(m)
         self.assert_shape_gamma(m)
+        self.assert_shape_vega(m)
+        self.assert_shape_theta(m)
         self.assert_shape_price(m)
         self.assert_shape_forward(m)
diff --git a/tests/nn/bs/test_bs.py b/tests/nn/modules/bs/test_bs.py
similarity index 100%
rename from tests/nn/bs/test_bs.py
rename to tests/nn/modules/bs/test_bs.py
diff --git a/tests/nn/bs/test_european.py b/tests/nn/modules/bs/test_european.py
similarity index 92%
rename from tests/nn/bs/test_european.py
rename to tests/nn/modules/bs/test_european.py
index 9ee25b50..5023d1ea 100644
--- a/tests/nn/bs/test_european.py
+++ b/tests/nn/modules/bs/test_european.py
@@ -181,9 +181,9 @@ def test_gamma_1(self):
 
     def test_gamma_2(self):
         m = BSEuropeanOption(call=False)
-        with pytest.raises(ValueError):
-            # not yet supported
-            _ = m.gamma(torch.tensor(0.0), torch.tensor(1.0), torch.tensor(0.2))
+        result = m.gamma(torch.tensor(0.0), torch.tensor(1.0), torch.tensor(0.2))
+        expect = torch.full_like(result, 1.9847627374)
+        assert_close(result, expect)
 
     def test_price_1(self):
         m = BSEuropeanOption()
@@ -217,6 +217,17 @@ def test_vega(self):
         expect = torch.tensor([0.1261, 0.1782, 0.2182])
         assert_close(result, expect, atol=1e-3, rtol=0)
 
+    def test_theta(self):
+        input = torch.tensor([[0.0, 0.1, 0.2], [0.0, 0.2, 0.2], [0.0, 0.3, 0.2]])
+        m = BSEuropeanOption(strike=100)
+        result = m.theta(
+            log_moneyness=input[..., 0],
+            time_to_maturity=input[..., 1],
+            volatility=input[..., 2],
+        )
+        expect = torch.tensor([-12.6094, -8.9117, -7.2727])
+        assert_close(result, expect, atol=1e-3, rtol=0)
+
     def test_example(self):
         from pfhedge.instruments import BrownianStock
         from pfhedge.instruments import EuropeanOption
@@ -235,5 +246,7 @@ def test_shape(self):
         m = BSEuropeanOption()
         self.assert_shape_delta(m)
         self.assert_shape_gamma(m)
+        self.assert_shape_vega(m)
+        self.assert_shape_theta(m)
         self.assert_shape_price(m)
         self.assert_shape_forward(m)
diff --git a/tests/nn/bs/test_european_binary.py b/tests/nn/modules/bs/test_european_binary.py
similarity index 82%
rename from tests/nn/bs/test_european_binary.py
rename to tests/nn/modules/bs/test_european_binary.py
index edf5cf0b..6440d202 100644
--- a/tests/nn/bs/test_european_binary.py
+++ b/tests/nn/modules/bs/test_european_binary.py
@@ -159,10 +159,17 @@ def test_forward(self):
 
     def test_delta(self):
         m = BSEuropeanBinaryOption()
-        result = m.delta(0.0, 0.1, 0.2)
+        result = m.delta(torch.tensor(0.0), torch.tensor(0.1), torch.tensor(0.2))
         expect = torch.tensor(6.3047)
         assert_close(result, expect, atol=1e-4, rtol=1e-4)
 
+        with pytest.raises(ValueError):
+            # not yet supported
+            m = BSEuropeanBinaryOption(call=False)
+            result = m.gamma(torch.tensor(0.0), torch.tensor(0.1), torch.tensor(0.2))
+            # expect = torch.tensor(-6.30468)
+            # assert_close(result, expect, atol=1e-4, rtol=1e-4)
+
     def test_gamma(self):
         m = BSEuropeanBinaryOption()
         result = m.gamma(torch.tensor(0.01), torch.tensor(1.0), torch.tensor(0.2))
@@ -173,6 +180,34 @@ def test_gamma(self):
             # not yet supported
             m = BSEuropeanBinaryOption(call=False)
             result = m.gamma(torch.tensor(0.0), torch.tensor(1.0), torch.tensor(0.2))
+            # expect = torch.tensor(-6.30468)
+            # assert_close(result, expect, atol=1e-4, rtol=1e-4)
+
+    def test_vega(self):
+        m = BSEuropeanBinaryOption()
+        result = m.vega(torch.tensor(0.0), torch.tensor(0.1), torch.tensor(0.2))
+        expect = torch.tensor(-0.06305)
+        assert_close(result, expect, atol=1e-4, rtol=1e-4)
+
+        with pytest.raises(ValueError):
+            # not yet supported
+            m = BSEuropeanBinaryOption(call=False)
+            result = m.vega(torch.tensor(0.0), torch.tensor(0.1), torch.tensor(0.2))
+            # expect = torch.tensor(0.06305)
+            # assert_close(result, expect, atol=1e-4, rtol=1e-4)
+
+    def test_theta(self):
+        m = BSEuropeanBinaryOption()
+        result = m.theta(torch.tensor(0.0), torch.tensor(0.1), torch.tensor(0.2))
+        expect = torch.tensor(0.0630)
+        assert_close(result, expect, atol=1e-4, rtol=1e-4)
+
+        with pytest.raises(ValueError):
+            # not yet supported
+            m = BSEuropeanBinaryOption(call=False)
+            result = m.theta(torch.tensor(0.0), torch.tensor(0.1), torch.tensor(0.2))
+            # expect = torch.tensor(-0.0630)
+            # assert_close(result, expect, atol=1e-4, rtol=1e-4)
 
     def test_price(self):
         m = BSEuropeanBinaryOption()
@@ -216,6 +251,8 @@ def test_shape(self):
 
         m = BSEuropeanBinaryOption()
         self.assert_shape_delta(m)
-        # self.assert_shape_gamma(m)
+        self.assert_shape_gamma(m)
+        self.assert_shape_vega(m)
+        self.assert_shape_theta(m)
         self.assert_shape_price(m)
         self.assert_shape_forward(m)
diff --git a/tests/nn/bs/test_lookback.py b/tests/nn/modules/bs/test_lookback.py
similarity index 76%
rename from tests/nn/bs/test_lookback.py
rename to tests/nn/modules/bs/test_lookback.py
index c9026fef..95c1bf5f 100644
--- a/tests/nn/bs/test_lookback.py
+++ b/tests/nn/modules/bs/test_lookback.py
@@ -12,6 +12,8 @@
 from ._utils import compute_delta
 from ._utils import compute_gamma
 from ._utils import compute_price
+from ._utils import compute_theta
+from ._utils import compute_vega
 
 
 class TestBSLookbackOption(_TestBSModule):
@@ -95,6 +97,68 @@ def test_check_gamma(self):
         expect = torch.tensor([0.0])
         assert_close(result, expect)
 
+    def test_check_vega(self):
+        m = BSLookbackOption()
+
+        # gamma = 0 for max --> +0
+        result = compute_vega(m, torch.tensor([[-10.0, -10.0, 0.1, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # gamma = 0 for max / spot --> +inf
+        result = compute_vega(m, torch.tensor([[0.0, 10.0, 0.1, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # delta = 0 for time --> +0
+        result = compute_vega(m, torch.tensor([[-0.01, -0.01, 1e-10, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        result = compute_vega(m, torch.tensor([[0.0, 0.01, 1e-10, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # delta = 0 for spot / k < 1 and volatility --> +0
+        result = compute_vega(m, torch.tensor([[-0.01, -0.01, 0.1, 1e-10]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        result = compute_vega(m, torch.tensor([[0.0, 0.01, 0.1, 1e-10]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+    def test_check_theta(self):
+        m = BSLookbackOption()
+
+        # gamma = 0 for max --> +0
+        result = compute_theta(m, torch.tensor([[-10.0, -10.0, 0.1, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # gamma = 0 for max / spot --> +inf
+        result = compute_theta(m, torch.tensor([[0.0, 10.0, 0.1, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # delta = 0 for time --> +0
+        result = compute_theta(m, torch.tensor([[-0.01, -0.01, 1e-10, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        result = compute_theta(m, torch.tensor([[0.0, 0.01, 1e-10, 0.2]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        # delta = 0 for spot / k < 1 and volatility --> +0
+        result = compute_theta(m, torch.tensor([[-0.01, -0.01, 0.1, 1e-10]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
+        result = compute_theta(m, torch.tensor([[0.0, 0.01, 0.1, 1e-10]]))
+        expect = torch.tensor([0.0])
+        assert_close(result, expect)
+
     def test_check_price(self):
         m = BSLookbackOption()
 
@@ -198,6 +262,8 @@ def test_shape(self):
         m = BSLookbackOption()
         self.assert_shape_delta(m)
         self.assert_shape_gamma(m)
+        self.assert_shape_vega(m)
+        self.assert_shape_theta(m)
         self.assert_shape_price(m)
         self.assert_shape_forward(m)
 
diff --git a/tests/test_autogreek.py b/tests/test_autogreek.py
index 8e5bfb71..7dada989 100644
--- a/tests/test_autogreek.py
+++ b/tests/test_autogreek.py
@@ -21,3 +21,15 @@ def pricer(volatility, coef):
     result = autogreek.vega(pricer, variance=variance, coef=coef)
     expect = 2 * coef * variance.sqrt()
     assert_close(result, expect)
+
+
+def test_theta():
+    def pricer(time_to_maturity, coef):
+        return coef * time_to_maturity.pow(2)
+
+    torch.manual_seed(42)
+    time_to_maturity = torch.randn(10).exp()  # make it positive
+    coef = torch.randn(10)
+    result = autogreek.theta(pricer, time_to_maturity=time_to_maturity, coef=coef)
+    expect = -2 * coef * time_to_maturity
+    assert_close(result, expect)