Optimal output synchronization of multi-agent leader-follower systems with unknown nonlinear dynamics is considered. The agents are assumed heterogeneous so that the dynamics may be nonidentical. A distributed observer is designed to estimate the leader state for each agent. A discounted performance function is defined for each agent, and an augmented Hamilton-Jacobi-Bellman (HJB) equation is derived to find its minimal value. The HJB solution depends on the trajectories of the local state and the distributed observer state. A control protocol based on the HJB solution assures that the synchronization error goes to zero locally asymptotically fast for all agents. The proposed approach has two main advantages compared to standard output synchronization methods. First, it is optimal in the sense that it not only makes the steady-state synchronization error zero, but also minimizes the transient error. Second, it does not require the explicit solution to the output regulator equations, though the HJB solutions implicitly provide optimal solutions to them. Finally, a reinforcement learning technique is used to learn the optimal control protocol for each agent without requiring any knowledge of the agents or the leader dynamics. Simulation studies on a notional multi-agent system validate the proposed approach.