I realized that it would be time consuming to post the final answer with all those supermathematics details here.

I would like to have a brief summary of the answer.

Let us denote the ($\mathcal{N}=1$) super Poincare group by $S\Pi$, and denote the set of real (complex) even Grassmann numbers by $\mathbb{R}_{c}$ ($\mathbb{C}_{c}$) and the set of real (complex) odd Grassmann numbers by $\mathbb{R}_{a}$ ($\mathbb{C}_{c}$). We define the superspace $\mathbb{R}^{4|4}$ as $(\mathbb{R}_{c})^{4}\times (\mathbb{R}_{a})^{4}$.

In addition, we define the group $SL(2,\mathbb{C}_{c})$ as the special linear group with entries being $\mathbb{C}_{c}$ numbers. Similarly, the group $SO(3,1|\mathbb{R}_{c})^{+}$ is the Lorentz group with entries being commutative real Grassmann numbers.

It can be shown that $SL(2,\mathbb{C}_{c})/\mathbb{Z}_{2}=SO(3,1|\mathbb{R}_{c})^{+}$. The full super Poincare group is $S\Pi=T^{4|0}(\mathbb{R}_{c})\cup SO(3,1|\mathbb{R}_{c})^{+}\cup T^{0|4}(\mathbb{C}_{a})$. i.e. $S\Pi$ has three subgroups which are identified as translation on even coordinates, Lorentz transformation on even and odd coordinates, and translation on odd coordinates, respectively.

It follows that we can prove that the superspace $\mathbb{R}^{4|4}$ is isomorphic the coset space $S\Pi/SO(3,1|\mathbb{R}_{c})^{+}$.

Remark: Since we are talking about $\mathcal{N}=1$ super Poincare group, we should use the identification $\mathbb{C}_{a}=\mathbb{R}_{a}\times\mathbb{R}_{a}$, and so $(\mathbb{R}_{c})^{4}\times(\mathbb{C}_{a})^{2}=\mathbb{R}^{4|4}$. In this way, it is possible to regard the odd coordinates as Majorana spinors, which will be clarified as follows.

An arbitrary group element of $S\Pi$ is of the form

$g(b,\epsilon,\bar{\epsilon},K)=\exp\left[-i(-b^{a}P_{a}+\frac{1}{2}K^{ab}j_{ab}+\epsilon^{\alpha}q_{\alpha}+\bar{\epsilon}^{\dot{\alpha}}\bar{q}_{\dot{\alpha}})\right]$,

where $\left\{P_{a},j_{ab},q_{\alpha},\bar{q}_{\dot{\alpha}}\right\}$ is a chosen basis of the Poincare superalgebra $\mathcal{G}$ over $\mathbb{C}$, the parameters $b^{a}$ and $K^{ab}$ are Grassmann even numbers, and $\epsilon_{\alpha}$ and $\bar{\epsilon}_{\dot{\alpha}}=(\epsilon_{\alpha})^{\ast}$ are Grassmann odd numbers.

The $\mathcal{N}=1$ Poincare superalgebra is a $\mathbb{Z}_{2}$-graded Lie algebra over $\mathbb{C}$, such that

$[P_{a},P_{b}]=0$,

$[j_{ab},P_{c}]=i\eta_{ac}P_{b}-i\eta_{bc}P_{a}$,

$[j_{ab},j_{cd}]=i\eta_{ac}j_{bd}-i\eta_{ad}j_{bc}+i\eta_{bd}j_{ac}-i\eta_{bc}j_{ad}$,

$[j_{ab},q_{\alpha}]=i(\sigma_{ab})^{\beta}_{\alpha}q_{\beta}$,

$[j_{ab},\bar{q}^{\dot{\alpha}}]=i(\bar{\sigma}_{ab})^{\dot{\alpha}}_{\dot{\beta}}\bar{q}^{\dot{\beta}}$,

$[P_{a},q_{\alpha}]=0$, $[P_{a},\bar{q}^{\dot{\alpha}}]=0$,

$\left\{q_{\alpha},q_{\beta}\right\}=0$,

$\left\{\bar{q}_{\dot{\alpha}},\bar{q}_{\dot{\beta}}\right\}=0$,

$\left\{q_{\alpha},\bar{q}_{\dot{\alpha}}\right\}=(\sigma^{a})_{\alpha\dot{\alpha}}P_{a}$.

This is simply a $\mathbb{Z}_{2}$-graded vector space over $\mathbb{C}$, denoted as $\mathcal{G}={}^{0}\mathcal{G}\oplus{}^{1}\mathcal{G}$, with a graded Lie bracket. The even generators $P_{a}$ and $j_{ab}$ form the ordinary Poincare algebra $\mathfrak{g}$ with ordinary Lie bracket $[\cdots,\cdots]$. The odd generators $q_{\alpha}$ and $\bar{q}_{\dot{\alpha}}$ are defined to anticommute with themselves and form a $\mathfrak{g}$-module via adjoint action.

While considering a generalization of ordinary Lie group to contain odd degrees of freedom, it is hopeless to use the exponential of the Poincare superalgebra since we do not have a generalization of Baker-Hausdorff formula for graded Lie brackets.

Instead, we may consider its Grassmann-shell, denoted as $\mathbf{G}$, which is defined as $\mathbf{G}=\Lambda_{\infty}\otimes\mathcal{G}$, where $\Lambda_{\infty}$ is the Grassmann algebra. This is a $\mathbb{Z}_{2}$-graded vector space over over $\Lambda_{\infty}$. i.e. $\mathbf{G}={}^{0}(\Lambda_{\infty}\otimes\mathcal{G})\oplus{}^{1}(\Lambda_{\infty}\otimes\mathcal{G})$. The first (second) summand is called the even (odd) part of the Grassmann-shell. In addition, we enforce a graded scalar multiplication rule such that any $\mathbb{C}_{a}$-valued Grassmann number anticommutes with odd vectors in the Grassmann-shell.

It is easy to show that the even part of the Grassmann-shell is closed under ordinatry Lie bracket multiplication $[\cdots,\cdots]$. Using the even part of the Grassmann-shell, ${}^{0}(\Lambda_{\infty}\otimes\mathcal{G})$, with its ordinary (antisymmetric) commutator $[\cdots,\cdots]$, we can obtain a generalization of ordinary Lie group, which is called the super Lie group, whose local coordinates are labelled with both commutative numbers and anticommutative numbers. Elements in the super Lie group can be obtained via the exponential of elements in the even part of the Grassmann-shell ${}^{0}(\Lambda_{\infty}\otimes\mathcal{G})$.

In our case, the even numbers $b_{a}$ and $K_{ab}$ are not necessarily real numbers, but can be any Grassmann-even numbers in $\mathbb{R}_{c}$ and odd numbers $\epsilon^{\alpha}$ and $\bar{\epsilon}^{\dot{\alpha}}$ are from $\mathbb{C}_{a}$. These Grassmann numbers play the role of local coordinates of the super Poincare group $S\Pi$.

Further, it can be shown that this superspace is an $S\Pi$-module, whose odd coordinates transform as Majorana spinors under the $SL(2,\mathbb{C}_{c})$ action.

The transformations are given as follows

$x^{\prime a}=(\exp{K})^{a}_{b}x^{b}$,

$\theta^{\prime}_{\alpha}=\left(\exp(\frac{1}{2}K^{ab}\sigma_{ab})\right)^{\beta}_{\alpha}\theta_{\beta}$,

$\bar{\theta}^{\prime\dot{\alpha}}=\left(\exp(\frac{1}{2}K^{ab}\bar{\sigma}_{ab})\right)^{\dot{\alpha}}_{\dot{\beta}}\bar{\theta}^{\dot{\beta}}$,

$x^{\prime a}=x^{a}+b^{a}$,

$\theta^{\prime\alpha}=\theta^{\alpha}+\epsilon^{\alpha}$,

$\bar{\theta}^{\prime\dot{\alpha}}=\bar{\theta}^{\dot{\alpha}}+\bar{\epsilon}^{\dot{\alpha}}$,

for arbitrary $\left\{x^{0},x^{1},x^{2},x^{3}|\theta^{\alpha},\bar{\theta}^{\dot{\alpha}}\right\}$ and $\left\{x^{\prime 0},x^{\prime 1},x^{\prime 2},x^{\prime 3}|\theta^{\prime\alpha},\bar{\theta}^{\prime\dot{\alpha}}\right\}\in\mathbb{R}^{4|4}$.

From the above expression, we see that we can indeed view the odd part of the superspace as Majorana spinors.

The physical intuition behind these super-generalizations is as follows.

1. Usually, people would like to treat the classical bosonic fields as complex-valued function. We may also treat the classical fermionic fields as complex-valued spinors. In most textbooks, we consider the unitary irreducible representation of the Poincare group as the 1-particle state. Knowing the transformation rule of 1-P states under action of the Poincare group, it is possible to recover the field equation satisfied by the classical fields. This is shown in the paper "Unitary Representations of the inhomogeneous Lorentz Group and their Significance in Quantum Physics" by Norbert Straumann.

http://arxiv.org/pdf/0809.4942v1.pdf

Following the procedure, we obtain complex-valued Dirac spinor fields. However, in quantum physics, taking the classical limit $\hbar\rightarrow 0$, the commutators (anticommutators)

$[\hat{\phi}^{I}(t,\vec{x}),\hat{\phi}^{J}(t,\vec{y})]_{\pm}=o(\hbar)$ leads to either Grassmann-even or Grassmann-odd classical fields defined as

$\phi^{I}(x)$: $\mathbb{R}^{4}\rightarrow\mathbb{C}^{p|q}$.

In this classical picture, the canonical quantization is much clearer: we replace commutative Grassmann numbers by elements of the Weyl algebra, and replace anti-commutative Grassmann numbers by elements of the Clifford algebra.

2. The picture of complex-valued classical fields is correct as long as we consider only non-interacting free fields. Let us consider the classical electrodynamics, given by the following action

$S=-\int d^{4}x\left\{F^{ab}F_{ab}+\bar{\Psi}{\not}D\Psi\right\}$,

the non-homogeneous equation of motion $\partial_{a}F^{ab}=\bar{\Psi}\gamma^{\mu}\Psi$ shows that the bosonic fields $A_{a}$ are propagated by the charges carried by fermionic 'bi-fields', while the homogeneous equation of motion $\partial_{[a}F_{bc]}=0$ tells us that the bosonic fields $A_{a}$ has free-dynamics. This inevitably forces us to re-consider about a better definition of classical fields.

To show this, let us review the Grassmann algebra $\Lambda_{\infty}$ generated by infinitely many generators $\xi^{1}$,$\cdots$, satisfying $\xi^{i}\xi^{j}+\xi^{j}\xi^{i}=0$. Any element $z\in\Lambda_{\infty}$ can be represented as

$z=z_{B}+z_{S}=z_{B}+\sum_{k=1}^{\infty}\sum_{i_{1}\cdots i_{k}}C_{i_{1}\cdots i_{k}}\xi^{i_{1}}\cdots\xi^{i_{k}}$,

where $z_{B}$, $C_{i_{1}\cdots i_{k}}\in\mathbb{C}$. We call $z_{B}$ the body of $z$ and $z_{S}$ the soul of $z$.

Inspired by the above definitions, it is much clearer if we think of the classical bosonic field $A_{a}$ as Grassmann-even-valued fields, whose body as a $\mathbb{R}$-valued function propagates in spacetime freely and soul is propagated by the charge carried by the Grassmann even source.

In the classical limit, one would also expect that all the spin degrees of freedom vanish as the planck constant $\hbar$ goes to $0$. The reason we treat classical fields as Grassmann-valued functions is to allow us to have a formal construction with a consistent classical action with interactions. To some extend, we should say that the electromagnetic radiation is purely a quantum effect, which is absent in the classical limit.

This construction also explains why we use Grassmann-odd-valued spinors for the path-integral formalism of fermions.

In this new picture of classical field theory, both classical bosonic fields and classical fermionic fields are unified as

$\phi^{I}:$ $\mathbb{R}_{c}^{4}\rightarrow\mathbb{C}^{p|q}$.

In this new picture, the classical Dirac fields $\bar{\Psi}$ and $\Psi$ should be independent elements of the set of Grassmann-valued fields. As a result, it is inappropriate to define the relation $\bar{\Psi}=\Psi^{\dagger}\gamma^{0}$. But still, we can have the involution of the complex Grassmann algebra defined as

$(zw)^{\ast}=(w)^{\ast}(z)^{\ast}$.

I am not sure how to find the rule of involution relating spinors and anti-spinors so that the classical Dirac Lagrangian is real.

In addition, the bonus benefit of this picture is that it allows us to have the 'classical' mass term of Majorana fermions.

All the details can be found from the book 'Ideas and Methods of Supersymmetry and Supergravity or A Walk Through Superspace' by Ioseph L Buchbinder and Sergei M Kuzenko.

Please leave a comment if you have any idea how to show that the Dirac Lagrangian takes values in $\mathbb{R}_{c}$.