In order to make reliable predictions with controlled uncertainties for a wide range of nuclear phenomena, a theoretical bottom-up approach, by which hadrons emerge from the underlying theory of strong interactions, quantum chromodynamics (QCD), is desired. The strongly interacting quarks and gluons at low energies are responsible for all the dynamics of nucleons and their clusters, the nuclei. The theoretical framework and the combination of analytical and numerical tools used to carry out a rigorous non-perturbative study of these systems from QCD is called lattice QCD. The result of a lattice QCD calculation corresponds to that of nature only in the limit when the volume of the spacetime is taken to infinity and the spacing between discretized points on the lattice is taken to zero. A better understanding of these discretization and volume effects, not only provides the connection to the infinite-volume continuum observables, but also leads to optimized calculations that can be performed with available computational resources. This thesis includes various formal developments in this direction, along with proposals for novel improvements, to be used in the upcoming LQCD studies of nuclear and hadronic systems. As the space(time) is discretized on a (hyper)cubic lattice in (most of) lattice QCD calculations, the lattice correlation functions are not fully rotationally invariant. This is known to lead to mixing between operators (those interpolating the states or inserting external currents) of higher dimensions with those of lower dimensions where the coefficients of latter diverge with powers of inverse lattice spacing, a, as the continuum limit is approached. This issue has long posed computational challenges in lattice spectroscopy of higher spin states, as well as in the lattice extractions of higher moments of hadron structure functions. We have shown, through analytical perturbative investigations of field theories, including QCD, on the lattice that a novel choice of operators, smeared over several lattice sites and deduced from a continuum angular momentum, has a smooth continuum limit. The scaling of the lower dimensional operators is proven to be no worse than a squared, explaining the success of recent numerical studies of excited state spectroscopy of hadrons with similar choices of operators. These results are presented in chapter 2 of this thesis. Due to Euclidean nature of lattice correlation function, the physical scattering parameters must be obtained via an analytical continuation to Minkowski spacetime. However, this continuation is practically impossible in the infinite-volume limit of lattice correlation function except at the kinematic thresholds. A formalism due to Luscher overcomes this issue by making the connection between the finite-volume spectrum of two interacting particles and their infinite-volume scattering phase shifts. We have extended the Luscher methodology, using an effective field theory approach, to the two-nucleon systems with arbitrary spin, parity and total momentum (in the limit of exact isospin symmetry) and have studied its application to the deuteron system, the lightest bound states of the nucleons, by careful analysis of the finite-volume symmetries. A proposal is presented that enables future precision lattice QCD extraction of the small D/S ratio of the deuteron that is known to be due to the action of non-central forces. By investigating another scenario, we show how significant volume improvement can be achieved in the masses of nucleons and in the binding energy of the deuteron with certain choices of boundary conditions in a lattice QCD calculation of these quantities. These results are discussed in chapters 3, 4 and 5. In order to account for electromagnetic effects in hadronic systems, lattice QCD calculations have started to include quantum electrodynamic (QED). These effects are particularly interesting in studies of mass splittings between charged and neutral members of isospin multiplets, e.g. neutral and charged pions. Due to the infinite range of QED interactions large volume effects plaque these studies. Using a non-relativistic effective theory for electromagnetic interactions of hadrons, we analytically calculate, and numerically estimate, the first few finite-volume corrections (up to 1 over L to the 4th power where L is the spatial extent of the volume) to the masses of hadrons and nuclei at leading order in the QED coupling constant, but to all orders in the short-distance strong interaction effects. These results are presented in chapter 6.